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Abstract: HERWIG is a general-purpose Monte Carlo event generator, which in- 
cludes the simulation of hard lepton-lepton, lepton-hadron and hadron-hadron scat- 
tering and soft hadron-hadron collisions in one package. It uses the parton-shower 
approach for initial- and final-state QCD radiation, including colour coherence effects 
and azimuthal correlations both within and between jets. This article updates the 
description of HERWIG published in 1992, emphasising the new features incorporated 
since then. These include, in particular, the matching of first-order matrix elements 
with parton showers, a more correct treatment of spin correlations and heavy quark 
decays, and a wide range of new processes, including many predicted by the Minimal 
Supersymmetric Standard Model, with the option of R-parity violation. At the same 
time we offer a brief review of the physics underlying HERWIG, together with details 
of the input and control parameters and the output data, to provide a self-contained 
guide for prospective users of the program. This version of the manual (version 3) 
is updated to HERWIG version 6.5, which is expected to be the last major release 
of Fortran HERWIG. Future developments will be implemented in a new C++ event 
generator, HERWIG++. 

Keywords: Supersymmetric Standard Model, QCD, Phenomenological Models^ 
Hadronic Collider^ . 
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1. Introduction 

HERWIG is a general-purpose event generator for high-energy processes, with par- 
ticular emphasis on the detailed simulation of QCD parton showers. The program 
provides a full simulation of hard lepton-lepton, lepton-hadron and hadron-hadron 
scattering and soft hadron-hadron collisions in a single package, and has the following 
special features: 

• Initial- and final-state QCD jet evolution with soft gluon interference taken 
into account via angular ordering; 

• Colour coherence of (initial and final) partons in all hard subprocesses, includ- 
ing the production and decay of heavy quarks and supersymmetric particles; 

• Azimuthal correlations within and between jets due to gluon interference and 
polarization; 

• A cluster model for jet hadronization based on non-perturbative gluon splitting, 
and a similar cluster model for soft and underlying hadronic events; 

• A space-time picture of event development, from parton showers to hadronic de- 
cays, with an optional colour rearrangement model based on space-time struc- 
ture. 

Several of these features were already present in HERWIG version 5.1 and were de- 
scribed accordingly in some detail in rcf. [1]. This was in turn based on the work of 
refs. [2-14]. In the present article we concentrate therefore on the new features in- 
corporated into HERWIG since 1992. The most important of these are the matching 
of first-order matrix elements with parton showers [15-20], a more correct treatment 
of spin correlations [21] and heavy quark decays, again including matrix element 
matching [19], and a wide range of new hard processes. In particular HERWIG now 
includes the production and decay of supersymmetric particles, with or without the 
assumption of R-parity conservation [22]. 

The precise HERWIG version described here is 6.500, which we shall generally 
refer to as "version 6" in the following. 

Let us recall briefly the main features of a generic hard, i.e. high momentum 
transfer, process of the type simulated by HERWIG. It can be divided notionally into 
four components, corresponding roughly to increasing scales of distance and time: 

1. Elementary hard subprocess. A pair of incoming beam particles or their con- 
stituents interact to produce one or more primary outgoing fundamental ob- 
jects. This can be computed exactly to lowest order in perturbation theory. 
The hard momentum transfer scale Q together with the colour flow of the 
subprocess set the boundary conditions for the initial- and final-state parton 
showers, if there are any. 
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2. Initial- and final-state parton showers. A parton constituent of an incident 
beam hadron with low spacelike virtuality can radiate timelike partons. In the 
process it decreases its energy to a fraction x of that of the beam and increases 
its spacelike virtuality, which is bounded in absolute value by the scale Q of 
the hard subprocess. This initial-state emission process leads to the evolution 
of the structure function F{x, Q) of the incident hadron. On a similar time- 
scale, an outgoing parton with large timelike virtuality can generate a shower 
of partons with lower virtuality. The amount of emission depends on the upper 
limit on the virtuality of the initiating parton, which is again controlled by the 
momentum transfer scale Q of the hard subprocess. Timelike partons from the 
initial-state emission may in turn initiate parton showering. The coherence of 
soft gluon emission from different parton showers is controlled by the colour 
flow of the subprocess. Within the showers, it is simulated by angular ordering 
of successive emissions. 

3. Heavy object decays. Massive produced objects such as top quarks, electroweak 
gauge and Higgs bosons, and possibly non-Standard Model (e.g. SUSY) parti- 
cles, can decay on time-scales that may be shorter than or comparable to that 
of the QCD parton showers. Depending on their nature and the decay mode, 
they may also initiate parton showers before and/or after decaying. 

4. Hadronization process. In order to construct a realistic simulation one needs 
to combine the partons into hadrons. This applies to the constituent partons 
of incoming hadronic beams as well as to the outgoing products of parton 
showering, which give rise to hadronic jets. This hadronization process takes 
place at a low momentum transfer scale, for which the strong coupling as 
is large and perturbation theory is not applicable. In the absence of a firm 
theoretical understanding of non-perturbative processes, it must be described 
by a phenomenological model. The model adopted in HERWIG is intended 
to disrupt as little as possible the event structure established in the parton 
showering phase. Showering is terminated at a low scale, Qq < 1 GeV, and 
the preconfinement property of perturbative QCD [23,24] is exploited to form 
colour-neutral clusters [4] which decay into the observed hadrons. Initial-state 
partons are incorporated into the incoming hadron beams through a soft, non- 
perturbative "forced branching" phase of spacehke showering. The remnants 
of incoming hadron, i.e. those constituent partons which do not participate in 
the hard subprocess, undergo a soft "underlying event" interaction modelled 
on soft minimum bias hadron-hadron collisions. 

After a brief, practical section on the use of HERWIG, in the following sections 
we discuss each of these components in turn, concentrating on the new features since 
version 5.1. 
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2. Using HERWIG 



The latest version of the program, together with all relevant information, is always 
available via the official HERWIG web page: 

http : //hepwww . rl . ac . uk/ theory/ seymour/herwig/ 

which is also mirrored at CERN: 

http : //home . cern . ch/ seymour/herwig/ 
The program is written in Fortran and the user has to modify the main program 



HWIGPR to generate the type and number of events required. See section for a 
sample main program. The program operates by setting up parameters in common 
blocks and then calling a sequence of subroutines to generate an event. Parameters 
not set explicitly in the block data HWUDAT or in HWIGPR are set to default values in the 
main initialisation routine HWIGIN. Output data are delivered in the LEP standard 
common block HEPEVT [25,26]. Note that all real variables accessible to the user, 
including those in HEPEVT, are of type DOUBLE PRECISION. 

Since version 6.3, to take account of the increased energy and complexity of inter- 
actions at LHC and future colliders, the default value of the parameter NMXHEP, which 
sets the array sizes in the standard /HEPEVT/ common block, has been increased to 
4000. 

To generate events the user must first set up the beam particle names PARTI, 
PART2 (type CHARACTER*8) and momenta PBEAMl, PBEAM2 (in GeV/c), a process code 
IPROC and the number of events required MAXEV. 

See section § for beams and processes available. 

All analysis of generated events (histogramming, etc.) should be performed by 
the user-provided routines HWABEG (to initialise the run), HWANAL (to analyse an event) 
and HWAEND (to terminate the run). 

A detailed event summary is printed out for the first MAXPR events (default 
MAXPR = 1). Setting IPRINT = 2 lists the particle identity codes, properties and 
decay schemes used in the program. 

The programming language is standard Fortran 77 as far as possible. However, 
the following may require modification for running on some computers 

• Most common blocks are inserted by INCLUDE 'HERWIG65 . INC statements — 
the file HERWIG65 . INC is part of the standard program package. 

• Many common blocks are initialised by BLOCK DATA HWUDAT. Although BLOCK 
DATA is standard Fortran 77, it can cause linkage problems for some systems. 

• Variables of type DOUBLE COMPLEX are used, which may be called COMPLEX* 16 
on some systems. 
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In the INCLUDE file, common blocks have been regrouped so that all commons that 
were present in version 6.1 are unchanged. All new variables since version 6.1 are in 
separate common blocks, which are frozen after each version. 

It is anticipated that version 6.5, with minor modifications in the series 6.5xx, 
will be the last Fortran version of HERWIG. Substantial physics improvements will be 
included in the new C++ event generator HERWIG++ [27], to be released in 2003. 
Watch for progress on the HERWIG++ web page: 

http : / /www . hep . phy . cam . ac . uk/ ~gieseke/Herwig++/ 



3. Physics simulated by HERWIG 

3.1 Elementary hard subprocess 

In HERWIG version 6 there is a fairly large library of QCD, electroweak and super- 
symmetric elementary subprocesses. These are listed and discussed in section ^ 

The hard subprocess plays an important role in defining the phase space of any 
associated initial- and final-state parton showers. As discussed in ref. [1] and refer- 
ences therein, the parton showers are branching processes in which the branchings 
are ordered in angle from a maximum to a minimum value determined by the cut- 
off Qq. The maximum value is determined by the elementary subprocess and is due 
to interference among soft gluons. The general result [6,7,28] is that the initial and 
final branchings are approximately confined within cones around the incoming and 
outgoing partons from the elementary subprocess. For the branching of parton i, 
the aperture of the cone is defined by the direction of the other parton j that is 
colour-connected to i. The relation between soft gluon interference and the colour 
connection structure of the elementary subprocess leads to many detectable effects 
in hadronic final states. For recent examples, see e.g. refs. [29,30]. 

For a general process there are various contributions with different colour connec- 
tions. The HERWIG library of elementary subprocesses includes the separate colour 
connection contributions. In general there is some ambiguity in the separation of 
contributions which are suppressed by inverse powers of the number of colours A^^^. 
In earlier versions of HERWIG, these sub-leading terms were divided amongst the 
various colour connection contributions according to the recipe in ref. [7]. In ver- 
sion 6 the following improved prescription [31] has been followed, for both the QCD 
and SUSY QCD subprocesses. The matrix element squared for each colour 

connection is computed in the limit Nc —>■ oo and the corresponding contribution is 
defined as 

where | A^p represents the sum over all colour connections including terms sub-leading 
in Nc- This ensures that each colour connection contribution is positive-definite and 
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has the correct pole structure, while the sum of contributions yields the exact (Born- 
level) cross section. 

Parton emission into phase space regions which are outside the above-mentioned 
angular cones, called in the following "dead zones", do not contribute to leading 
order and often not even to next-to-leading order. However, for a more complete 
description of the event it is also necessary to take into account emission into these 
"dead zones" , which may be done using exact fixed-order matrix elements (see the 
following sections). 

Another important function of the elementary subprocess is to set up the po- 
larizations of any clcctroweak bosons, heavy (e.g. SUSY) particles or gluon jets that 
may be involved. These polarizations give rise to angular asymmetries and corre- 
lations in boson decays and jet fragmentation. They are included in HERWIG for 
many of the subprocesses provided, using the approach of refs. [9, 10,21] to generate 
all correlations in decays, and in jet fragmentation to leading-logarithmic accuracy. 

3.2 Parton showers 
3.2.1 Final-state showers 

Final-state parton showering in HERWIG is generated by a coherent branching algo- 
rithm with the following properties: 

1. The energy fractions are distributed according to the leading-order Dokshitzer- 
Gribov-Lipatov-Altarelli-Parisi (DGLAP) splitting functions. 

2. The full available phase space is restricted to an angular-ordered region. Such 
a restriction is the result of interference and correctly takes important infrared 
singularities into account. At each branching, the angle between the two emit- 
ted partons is smaller than that of the previous branching. 

3. The emission angles arc distributed according to the Sudakov form factors, 
which sum the virtual corrections and unresolved real emissions. The Sudakov 
form factor normalizes the branching distributions to give the probabilistic 
interpretation needed for a Monte Carlo simulation. This fact is a consequence 
of unitarity and of the infrared finiteness of inclusive quantities. 

4. The azimuthal angular distribution in each branching is determined by two 
effects: 

(a) for a soft emitted gluon the azimuth is distributed according to the eikonal 
dipole distribution [7]; 

(6) for non-soft emission one finds azimuthal correlations due to spin effects. 
See [9, 10] for the method used to implement these correlations in full, to 
leading coUinear logarithmic accuracy, in HERWIG. 
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5. In each branching the scale of as is the relative transverse momentum of the 
two daughter partons. 



6. In the case of heavy flavour production the mass of the quark modifies the 
angular-ordered phase space. The most important effect is that the soft radi- 
ation in the direction of the heavy quark is depleted. One finds that emission 
within an angle of order M/E is suppressed, M and E being the mass and 
energy of the heavy quark: this is known as the "dead cone" [32]. 

Specifically, the HERWIG parton shower evolution is done in terms of the parton 
energy fraction z and an angular variable ^. In the parton splitting i — > jfc, Zj = 
Ej/Ei and = {pj ■ pk)/{EjEk). Thus C,jk — for massless partons at small 
angles. The values of z are chosen according to the DGLAP splitting functions and 
the distribution of ^ values is determined by the Sudakov form factors. Angular 
ordering implies that each ^ value must be smaller than the ^ value for the previous 
branching of the parent parton. 

The parton showers are terminated as follows. For partons there is a cutoff of the 
form Qi = rrii + Qq, where mj (z = 1, . . . , 6 for d, u, s, c, b, t) is set by the relevant 
mass parameter RMASS(i) and Qq is set by the quark and gluon virtuality cutoff 
parameters VQCUT and VGCUT (see section ^). Showering from any parton stops when 
a value of ^ below Q]/Ef is selected for the next branching. For heavy quarks, the 
condition ^ > Qf/Ef corresponds to the "dead cone" mentioned above. At this point 
the parton is put on mass-shell or given a small non-zero effective mass in the case 
of gluons.^ Working backwards from these on-shell partons, one can now construct 
the virtual masses of all the internal lines of the shower and the overall jet mass, 
from the energies and opening angles of the branchings. Finally one can assign the 
azimuthal angles of the branchings, including EPR-type correlations (from Einstein- 
Podolski- Rosen [33]), and deduce completely all the 4- momenta in the shower. 

Next the parton showers are used to replace the (on mass-shell) partons that 
were generated in the original hard process. This is done in such a way that the 
jet 3-momenta have the same directions as the original partons in the cm. frame of 
the hard process, but they are boosted to conserve 4-momentum taking into account 
their extra masses. 

The main improvements in the final-state emission algorithm of HERWIG ver- 
sion 6, relative to version 5.1, are as follows. 

The Sudakov form factors can be calculated using the one-loop or two- loop as, 
according to the variable SUDORD (default = 1). The parton showering still incor- 
porates the two-loop as in either case but if SUDORD = 1 this is done using a veto 
algorithm, whereas if SUDORD = 2 no vetoes are used in the final-state evolution. The 



usefulness of this option is discussed briefly in section S.2 



^The quark mass parameters should also be thought of as effective or constituent masses rather 
than current quark masses. 
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Matrix element corrections have been introduced into final-state parton showers 
in e^e~ and deep inelastic processes [15,16], in heavy flavour decays [19] and in 
Drell-Yan processes [20] (see section p.2.3|) . 

In HERWIG, the angular-ordering constraint, which is derived for soft gluon 
emission, is applied to all parton shower vertices, including g qq. In versions before 
6.1, this resulted in a severe suppression (an absence in fact) of configurations in 
which the gluon energy is very unevenly shared between the quarks. For light quarks 
this is irrelevant, because in this region one is dominated by gluon emissions, which 
are correctly treated. However, for heavy quarks this energy sharing (or equivalently 
the quarks' angular distribution in their rest frame) is a directly measurable quantity 
and was badly described. Related to this was an inconsistency in the calculation of 
the Sudakov form factor for g — > qq. This was calculated using the entire allowed 
kinematic range (with massless kinematics) for the energy fraction, < z < 1, while 
the z distribution generated was actually confined to the angular-ordered region, 
z,l-z > m/E^. 

In version 6, these defects are corrected as follows. We generate the E,C, and z 
values for the shower as before. We then apply an a posteriori adjustment to the 
kinematics of the g ^ qq vertex during the kinematic reconstruction. At this stage, 
the masses of the q and q showers are known. We can therefore guarantee to stay 
within the kinematically allowed region. In fact, the adjustment we perform is purely 
of the angular distribution of the q and q showers in the g rest frame, preserving all the 
masses and the gluon four-momentum. Therefore we do not disturb the kinematics 
of the rest of the shower at all. 

Although this cures the inconsistency above, it actually introduces a new one: 
the upper limit for subsequent emission is calculated from the generated E, ^ and 
X values, rather than from the finally-used kinematics. This correlation is of NNLL 
importance, so we can formally neglect it. It would be manifested as an incorrect 
correlation between the masses and directions of the produced q and q jets. This 
is, in principle, physically measurable, but it seems less important than getting the 
angular distribution itself right. In fact the solution we propose maps the old an- 
gular distribution smoothly onto the new, so the sign of the correlation will still be 
preserved, even if the magnitude is wrong. Even with this modification, the HER- 
WIG kinematic reconstruction can only cope with particles that are emitted into the 
forward hemisphere in the showering frame. Thus one cannot populate the whole of 
kinematically-allowed phase space. Nevertheless, we find that this is usually a rather 
weak condition and that most of phase space is actually populated. 

Using this procedure, we find that the predicted angular distribution for sec- 
ondary b quarks at LEP energies is well-behaved, i.e. it looks reasonably similar to 
the leading-order result (1 + cos^^*), and has relatively small hadronization correc- 
tions. 

Real photon emission is included in timelike parton showering. The infrared 
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photon cutoff is VPCUT, which defaults to 0.4 GeV. Agreement with LEP data is sat- 
isfactory if showering is used together with the matrix element correction to produce 
photons in the back-to-back region. The results are insensitive to VPCUT variations in 
the range 0.1-1.0 GeV. Setting VPCUT greater than the total cm. energy switches off 
such emission. As an expedient way of improving the efficiency of photon final-state 
radiation studies, the electromagnetic coupling ALPHEM can be multiplied by a factor 
ALPFAC (default = 1) for all quark-photon vertices in jets, and in the 'dead zone' 
in e+e". Results at small photon-jet separation become sensitive to ALPFAC above 
about 5. 

3.2.2 Initial-state showers 

The theoretical analysis of initial-state showering is more complex than the final- 
state case. The most relevant parameters of the hard subprocess are the hard scale 
Q and the energy fraction x of the incoming parton after the emission of initial state 
radiation. For lepton-hadron processes x corresponds to the Bjorken variable, while 
for hadron-hadron processes x is related to Q^/s where s is the cm. energy squared. 

The main result is that for any value of x, even for x small [34], the initial-state 
emission process factorizes and can be described as a coherent branching process 
suitable for Monte Carlo simulations. The properties which characterize this process 
include all those discussed above for the final-state emission. However, in the initial- 
state case the angular-ordering restriction on the phase space applies to the angles 
9i between the directions of the incoming hadron and the emitted partons i. 

For large x, the coherent branching algorithm sums correctly [13] not only the 
leading but also the next-to-leading contributions. This accuracy allows us to iden- 
tify the relation between the QCD scale used in the Monte Carlo program and the 
fundamental parameter Aj^. This is achieved by using the one-loop Altarelli-Parisi 
splitting functions and the two-loop expression for as with the following universal 
relation between the scale parameter Aphys [13] used in the simulation and Aj^ (here, 
Nf is the number of flavours) 



Therefore a Monte Carlo simulation with next-to-leading accuracy can be used to 
determine A^^ from semi-inclusive data at large momentum fractions.^ 

In the case of small values of x, the initial state branching process has additional 
properties, which are not yet included fully in HERWIG. This was discussed in ref. [1] 
and the situation remains unchanged since version 5.1. 

The initial-state branching algorithm in HERWIG is of the backward evolution 
type. It proceeds from the elementary subprocess, at a hard scale set by colour co- 

^This applies also to final-state emission, i.e. to jet fragmentation at large values of the jet 
momentum fraction. 
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herence (see section |3.1| ), back to the hadron scale, set here by the spacehke cutoff 
parameter QSPAC. At this point there is a forced non-perturbative stage of branch- 
ing which ensures that the emitting parton fits smoothly with the valence parton 
distributions of the incoming hadron. 

Matrix-element corrections have been introduced into initial-state parton showers 
in deep inelastic [16] and Drell-Yan processes [20], as discussed in the following 
subsection. 

To avoid double-counting of hard parton emission, all radiation at transverse 
momenta greater than the hard process scale EMSCA is vetoed. In the case of initial- 
state radiation, this affects all events, while for final-state radiation it only affects 
those events in which the two jets have a rapidity difference of more than about 3.4. 

In the backward evolution of initial-state radiation for photons the "anomalous" 
branching gg ^ 7 is included. Variables AN0MSC(1,IBEAM) and AN0MSC(2,IBEAM) 
record the evolution scale and transverse momentum, respectively, at which an 
anomalous splitting was generated in the backward evolution of beam IBEAM. If zero, 
then no such splitting was generated. 

The treatment of forced branching of gluons and sea (anti-) quarks in backward 
evolution has been improved, by allowing it to occur at a random scale between the 
spacelike cutoff QSPAC and the infrared cutoff, instead of exactly at QSPAC as before. 
A new option ISPAC = 2 allows the freezing of structure functions at the scale QSPAC, 
while evolution continues to the infrared cutoff. The default, ISPAC = is equivalent 
to previous versions, in which perturbative evolution stops at QSPAC. 

The width of the (gaussian) intrinsic transverse momentum distribution of va- 
lence partons in incoming hadrons at scale QSPAC is set by the parameter PTRMS 
(default value zero). The intrinsic transverse momentum is chosen before the ini- 
tial state cascade is performed and is held fixed even if the generated cascade is 
rejected. This is done to avoid correlation between the amount of perturbative and 
non-perturbative transverse momentum generated. 

It is possible to completely switch off initial-state emission, by setting NOSPAC = 
. TRUE . , in which case only the forced splitting of non- valence partons is generated. 

3.2.3 Matrix-element corrections 

One of the new features of HERWIG 6 is the matching of first-order matrix elements 
with parton showers. 

The HERWIG parton showers are performed in the soft or collinear approximation 
and emission is allowed only in regions of the phase space satisfying the condition ^ < 
1 or ^ < for the final- (timelike branching) and initial-state (spacelike branching) 
radiation respectively, where ^ and z are the showering variables defined above. 

The emission is entirely suppressed inside the so-called dead zones (^ > 1 or 
> z^), corresponding to hard and/or large-angle parton radiation. According to 
the exact matrix elements, the radiation in the dead zones is suppressed, since it 
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is not soft or collinear logarithmically enhanced, but it is not completely absent as 
happens in the HERWIG standard shower algorithm. The HERWIG parton cascades 
need to be supplemented by matrix-element corrections for a full description of the 
physical phase space. 

The method of matrix-element corrections to the HERWIG parton showers is dis- 
cussed in [16,17]. The radiation in the dead zones is generated according to the exact 
first-order matrix element ('hard correction'); the shower in the already-populated 
region of the phase space is corrected by the use of the exact 0{as) amplitude any 
time an emission is capable of being the 'hardest so far' ('soft correction'^). 

By 'hardcst-so-far', we mean the radiation of a parton whose transverse momen- 
tum relative to the splitting one is larger than all those previously emitted. This is 
not always the first emission, as angular ordering does not necessarily imply ordering 
in transverse momentum. As shown in [16], if we corrected only the first emission, we 
would have problems in the implementation of the Sudakov form factor whenever a 
subsequent harder emission occurs, as we would find that the probability of hard 
radiation would depend on the infrared cutoff, which is clearly unphysical. Using the 
0{as) result for the hardcst-so-far emission in the already-filled phase space as well 
as in the dead zone allows one to have matching over the boundary of the dead zone 
itself. 

Since the fraction of events which receive a hard correction is typically small, we 
neglect multiple hard emissions in the dead zones and rely on the first-order result 
plus showering in those regions. 

Our method is quite different from the one used to implement matrix-element cor- 
rections in JETSET [35], where the parton shower probability is applied over the whole 
phase space and the first-order amplitude is used only to correct the first emission. 

Following these general prescriptions, matrix-element corrections have been im- 
plemented in some e~^e~ processes [15] (including e~^e~ WW/ZZ), deep inelastic 
Icpton scattering [18], top quark decay [19], and Drell-Yan processes [20]. Those for 
the process gg Higgs are now in progress [36] . 

The variables HARDME (default = .TRUE.) and SOFTME (default = .TRUE.) al- 
low respectively the application of hard and soft matrix-element corrections to the 
HERWIG parton cascades. 

3.3 Heavy flavour production and decay 

Heavy quark decays are treated as secondary hard subprocesses. Top quarks and 
any hypothetical heavier quarks always decay before hadronization. Heavy-flavoured 
hadrons are split into collinear heavy quark and spectator and the former decays in- 
dependently. After decay, parton showers may be generated from coloured decay 

^We point out that in the expression 'soft correction', 'soft' refers to the phase space where such 
corrections are applied and not to the amplitude, since we still use the 'hard' exact matrix element 
for the soft correction as well. 
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products, in the usual way. See ref. [11] for details of the treatment of colour coher- 
ence in these showers. 

In HERWIG version 6 matrix-element corrections to the simulation of top quark 
decays are available. The routine HWBTOP implements the hard corrections; HWBRAN 
has been modified to implement the soft corrections. Since the dead zone includes 
part of the soft singularity, a cutoff is required: only gluons with energy above GCUTME 
(default value 2 GcV) in the top rest frame are corrected. Physical quantities are not 
strongly dependent on GCUTME in the range 1 to 5 GeV, after the typical experimental 
cuts are applied. For more details sec ref. [19]. 

The structure of the program has been altered so that the secondary hard sub- 
process and subsequent fragmentation associated with each partonic heavy hadron 
decay appear separately in the event record. Thus top quark decays are treated 
individually as are any subsequent bottom hadron partonic decays. Note that the 
statement CALL HWDHOB, which deals with the decays of all heavy objects (including 
SUSY particles), must appear in the main program between the calls to HWBGEN and 
HWCFOR, in order to carry out any decays before hadronization. 

The partonic decay fractions of heavy quarks are specified in the decay tables 
like the decay modes of other particles. This permits different decays to be given 
to individual heavy hadrons. Changes to the decay table entries can be made on an 
event by event basis if desired. Partonic decays of charm hadrons and quarkonium 
states are also now supported. The order of the products in a partonic decay mode 
is significant. For example, if the decay is Q — > W + q ^ {f + f') + Q occurring inside 
a Qs hadron, the required orderings are: 

Q + s^{f + f') + {q + s) 

or {q + f) + if + s) ('colour rearranged') . 

In both cases the {V — A)"^ matrix element-squared is proportional to {pg •p2)(pi -ps), 
where pi etc. correspond to the ordering given. Decays of heavy-flavoured hadrons 
to exclusive non-partonic final states are also supported. No check is made against 
double counting from partonic modes. However this is not expected to be a major 
problem for the semi-leptonic and two-body hadronic modes supplied. 

The default masses of the c and b quarks have been lowered to 1.55 and 4.95 re- 
spectively: this corresponds to the mass of the lightest meson minus the u or d quark 
mass. This increases the number of heavy mesons, and hence total multiplicities, 
and slightly softens their momentum spectrum. The rate of photoproduced charm 
states increases and B-tt momentum correlations become smoother. The default top 
quark mass is 174.3 GeV/c^. The same value is used in the production and decay 
matrix elements and for all kinematics. Note that higher-order corrections are not 
fully included, and so the HERWIG top mass does not necessarily correspond to that 
defined in any particular rigorous scheme (e.g. the pole mass or the MS running 
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mass). However, since it is probably the decay kinematics that are most sensitive to 
this parameter, it should be close to the pole mass. See subsection [4.2. 1| for notes on 
the treatment of quark masses in various processes. 



3.4 Gauge and Higgs boson decays 

The total decay widths of the elec- 
troweak gauge bosons V = W, Z are 
specified by the input parameters GAMW 
and GAMZ. Their branching fractions to 
various final states are computed auto- 
matically from the other SM input pa- 
rameters. Which decays actually occur 
is controlled as follows. The variable 
MODBOS(i) controls the decay of the ith. 
gauge boson per event (table |I|). 

All entries of MODBOS default to 0. 
Bosons which are produced in pairs (i.e. Table 1: W, Z decays, 

from VV pair production, or Higgs de- 
cay) are symmetrized in MODBOS(i) and MODBOS(i + 1). For processes which directly 
produce gauge bosons, the event weight includes the branching fraction to the re- 
quested decay, but this is only true for Higgs production if decay to W~^W~ /Z^Z^ 
is forced (IPROC = 310, 311 but not 399, etc.). Users can thus force Z bb decays, 
with M0DB0S(2) = 7. For example, IPROC = 250, MODBOS(l) = 7, M0DB0S(2) = gives 
Z^Z^ production with one Z° decaying to bb. 

The spin correlations in the decays are handled in one of two ways: 

1. The diagonal members of the spin density matrix are stored in RHOHEP(i, IHEP), 
where i = 1,2,3 for helicity= i — 2 in the centre-of-mass frame of their produc- 
tion, for processes where this matrix is diagonal (i.e. there is no interference 
between spin states). 

2. The correlations in the decay are handled directly by the production routine 
where (1) is not possible. 

The processing of the parton showers in hadronic W and Z decays is handled 
in the rest frame of the vector boson if WZRFR is .TRUE, (the default), otherwise in 
the lab frame. In the latter case, which was the default in earlier versions, the initial 
cone angles of the showers depend on the velocity of the boson, which leads to a 
slight Lorentz non-invariance of decay distributions. 

The total decay width of the SM Higgs boson is computed from its input mass 
RMASS(201) and stored as GAMH. Its decay branching fractions are also computed and 
stored in BRHIG(I): I = 1-6 for dd,. ..,«;! = 7-9 for e+e",. . . , r+r"; I = 10,11,12 
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for W~^W^ , Z^Z^, 77. Non-SM Higgs bosons, on the other hand, such as those 
in super symmetric models, have to have their widths and decay tables provided as 
input data (see section |3.5.1| ). To avoid any ambiguity, the SM Higgs boson has a 
distinct identity code in HERWIG and is represented by the special symbol -ffgM- 

There are two choices for the treatment of the SM Higgs width, both controlled 
by the variable lOPHIG: 

lOPHIG = 2/ + J, 

where / and J are both zero or one. Whenever a Higgs boson is generated, its 
mass is chosen from a distribution that, for heavy SM Higgs bosons, can be rather 
broad. The choice of I makes a significant difference to the physical meaning of 
the distribution generated: for / = 0, the cross section corresponds to the tree level 
process containing an s-channel Breit-Wigner resonance for the Higgs boson with a 
running Higgs width. As discussed in [37], this neglects important contributions from 
interference with non-resonant diagrams and can violate unitarity at high energy, so 
1 = 1 (the default) uses the improved prescription of [37]. This replaces the s-channel 
propagator by an effective propagator that sums the interference terms to all orders. 
This increases the cross section below resonance and decreases it above, causing an 
overall increase in cross section. More details can be found in [37]. 

The variable J is a more technical parameter that does not affect the physical 
results, only the method by which they are generated: J = 1 (the default) generates 
the mass according to a fixed- width Breit-Wigner resonance, while J = biases the 
distribution more towards higher masses. In either case, the appropriate jacobian 
factor is included in the event weight, so that the physical cross section is independent 
of J. 

In all the above cases, the SM Higgs mass distribution is restricted to the range 
[mn - GAMMAX x rH,mH + GAMMAX x Th]. GAMMAX defaults to 10, but in the non- 
perturbative region mn ^ 1 TeV should be reduced to protect against poor weight 
distributions. These considerations do not affect the distribution noticeably for 
mn ^ 500 GeV, and GAMMAX can safely be increased if necessary. 

For a SUSY Higgs, the width is never large enough for unitarity to be violated and 
these issues are unimportant. In this case, the mass distribution is chosen according 
to a fixed-width Breit-Wigner resonance, like that of any other SUSY particle. 

The SM Higgs decays that can occur are normally controlled by the process code 
IPROC, as in IPROC = 300 + ID for example: ID= 1-6 for quarks, 7-9 for leptons, 
10/11 for W~^W~ /Z^Z^ pairs, and 12 for photons. In addition ID= gives quarks of 
all flavours, and ID= 99 gives all decays. For each process, the average event weight 
is the cross section in nb times the branching fraction to the requested decay. The 
branching ratios to quarks use the next-to-leading logarithm corrections, those to 
W^W~ /Z^Z^ pairs allow for one or both bosons being off mass-shell. 

All Higgs vertices include an optional enhancement factor to account for non-SM 
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and non-MSSM couplings. The amplitudes for all Higgs vertices are multiplied by 
the factor ENHANC(ID) where ID is the same as in IPROC = 300 + ID except the 77if 
'vertex' which is calculated from ENHANC(6) and ENHANC(IO) for the top and 
loops. This allows the simulation of the production of any chargeless scalar Higgs-like 
particle. Note however that pseudoscalar and charged Higgs bosons, and processes 
involving more than one Higgs particle (e.g. the decay —>■ h^Z) are not included 
this way (see section 

The array ENHANC(ID) is initialised as usual in HWIGIN. Note, however, that it 
will be overwritten if MSSM Higgs production is required by IPROC. In that case, as 
mentioned earlier, the Higgs widths and decay modes are simply read from an input 



particle data file (see section p.5.1 
3.5 Supersymmetry 



Particle 




Spin 


Particle 




Spin 


quark 




1/2 


squarks 


(lL,R 





charged lepton 


£ 


1/2 


charged sleptons 


h,R 





neutrino 


V 


1/2 


sneutrino 


V 





gluon 


9 


1 


gluino 


9 


1/2 


photon 


7 


1 


photino 


7 


1/2 


neutral gauge boson 




1 


zino 


Z 


1/2 


neutral Higgs bosons 







neutral Higgsinos 


-"1,2 


1/2 


charged gauge boson 




1 


wino 




1/2 


charged Higgs boson 







charged Higgsino 




1/2 


graviton 


G 


2 


gravitino 


G 


3/2 


W^, mix 


to form 2 chargino 


mass eigenstates 







7, Z, if°2 i^ix to form 4 neutralino mass eigenstates Xi;X25X3;X4 
t^, tji (and similarly b, r) mix to form the mass eigenstates ti, t2 

Table 2: SUSY particle content. 



HERWIG now includes the production and decay of superparticles, as given by the 
Minimal Supersymmetric Standard Model (MSSM) [22]. The mass spectrum and de- 
cay modes, being read from input files (see below), are completely general. The par- 
ticle content, listed in the following table (table ||), includes the gravitino/goldstino. 
For sparticles that mix, the subscripts label the mass eigenstates in the ascending 
order of mass. The two Higgs Doublet Model (2HDM) Higgs sector, intrinsic to 
the MSSM, is also included. The three neutral Higgs bosons are denoted by 
and y4°. 

The SUSY particle names in HERWIG are as shown in table ^. IDHW is the 
HERWIG identity code and IDPDG is the corresponding Particle Data Group code [38]. 
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IDHW 




NAME 


IDPDG 


IDHW 




NAME 


IDPDG 


401 




'SSDL ' 


1000001 


413 


dR 


'SSDR ' 


2000001 


402 


Ul 


'SSUL ' 


1000002 


414 


Ur 


'SSUR ' 


2000002 


403 


sl 


'SSSL ' 


1000003 


415 


sr 


'SSSR ' 


2000003 


404 


cl 


'SSCL ' 


1000004 


416 


cr 


'SSCR ' 


2000004 


405 


hi 


'SSBl ' 


1000005 


417 


b2 


'SSB2 ' 


2000005 


406 


ti 


'SSTl ' 


1000006 


418 


t2 


'SST2 ' 


2000006 



425 


Cl 


'SSEL- ' 


1000011 


437 


eR 


'SSER- ' 


2000011 


426 




'SSNUEL ' 


1000012 










427 




'SSMUL- ' 


1000013 


439 


V-R 


'SSMUR- ' 


2000013 


428 




'SSNUMUL ' 


1000014 










429 


n 


'SSTAUl- ' 


1000015 


441 




'SSTAU2- ' 


2000015 


430 




'SSNUTL ' 


1000016 











449 


9 


'GLUINO ' 


1000021 


458 


G 


'GRAVTINO' 


1000039 


450 




'NTLINOl ' 


1000022 


451 


xl 


'NTLIN02 ' 


1000023 


452 


xl 


'NTLIN03 ' 


1000025 


453 


xl 


'NTLIN04 ' 


1000035 


454 


xt 


'CHGIN01+' 


1000024 


455 


xt 


'CHGIN02+' 


1000037 


203 




'HIGGSLO ' 


26 


204 




'HIGGSHO ' 


35 


205 


AO 


'HIGGSAO ' 


36 


206 




'HIGGS+ ' 


37 



Table 3: SUSY particle names. 



Antiparticles generally appear in sequence after the corresponding particles, e.g. 
antisquarks dl-tl at IDHW = 407-412, d*j^-t*2 at 419-424. They have 'BR' added to 
the name, e.g. ' SSDLBR ' , or opposite charge, and negative PDG codes. A full list 
can be obtained using the print option IPRINT = 2 (see section 

Note that the HERWIG particle labelling of the lightest MSSM Higgs boson 
departs from the PDG recommendation: it is given PDG code 26, to avoid confusion 
with the SM Higgs boson (PDG code 25) in our implementation (specifically, in our 
use of the array ENHANC for the MSSM processes: see the relevant Higgs sections for 
more details). 

HERWIG does not contain any built-in models for SUSY scenarios beyond the 
MSSM, such as, Supergravity (SUGRA) or Gauge Mediated Symmetry Breaking 
(GMSB). In all cases the SUSY particle spectrum and decay tables must be provided 
just like those for any other particles. The subroutine HWISSP, if called, reads these 
from an input file. The production subprocesses are then generated by HWHESP, in 
lepton-antilepton collisions, HWHSSP, in hadron-hadron collisions, or by one of the 
fip production routines. The decays of the sparticles produced, as well as any top 
quarks or Higgs bosons, are then performed by HWDHOB. 
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3.5.1 Data input 



A package ISAWIG, see section dA, has been created to work with ISAJET [39] to 
produce a file containing the SUSY particle masses, hfetimes and decay modes. This 
package takes the outputs of the ISAJET SUGRA or general MSSM programs and 
produces a data file in a format that can be read into HERWIG by the subroutine 
HWISSP. In principle the user can produce a similar file provided that the correct 
format is used, as explained below. 

For the mixing terms of the MSSM lagrangian we follow the Haber-Kane [40,41] 
conventions, so that we differ from ISAJET on the sign for gaugino masses, the 
ordering and signs of the gaugino current eigenstates, the interchange of the rows 
and columns of the gaugino mixing matrices, and the sign of the neutral Higgs mixing 
angle a. 

In addition to the decay modes included in the ISAJET package ISAWIG allows 
for the possibility of violating R-parity and includes the calculation of all 2-body 
squark and slepton, and 3-body gaugino and gluino R-parity violating (fip) decay 
modes. 

It can happen that some of the SUSY particle decay modes generated by ISAJET 
are found to be kinematically forbidden in HERWIG, owing to the slightly different 
values assumed for the light quark masses. In this case a warning message is printed 
by HERWIG and these modes are deleted, the other branching ratios being rescaled 
accordingly. Such modes normally have negligible ISAJET branching ratios anyway, 
because of their tiny phase space. 

The input file organisation expected by HWISSP is as follows. First the SUSY 
particle and top quark masses and lifetimes (in seconds) are given according to their 
HERWIG identity codes IDHW, for example: 

65 

401 927.3980 0.74510E-25 

402 925.3307 0.74009E-25 
.... etc . 

That is, 

NSUSY = Number of SUSY + top particles 
IDHW, RMASS(IDHW), RLTIM(IDHW) 
repeated NSUSY times. 

Next each particle's decay modes together with their branching ratios and matrix 
element codes are given as, for example: 

6 

401 0.18842796E-01 450 1 
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401 0.32755006E-02 457 2 
6 

402 0.94147678E-02 450 2 
.... etc . 

That is, 

Number of decay modes for a given particle IDK 
IDK(IM), BRFRAC(IM), NME(IM), IDKPRD(l-5, IM) 
repeated for each mode IM 
all repeated for each particle (NSUSY times). 

The order in which the decay products appear is important in order to obtain ap- 
propriate showering and hadronization. The correct orderings are indicated in the 
table below (table g). 

New matrix element codes have been added for SUSY and Higgs decays: 

• NME=200, describing the 1^3 body heavy-quark decays via a virtual 
boson. A new function, HWDHWT, was introduced to this end. This can be used 
to emulate e.g. t — >• hH^{-^ ff) decays. 

• NME = 300 for three-body fip gaugino and gluino decays. 

The indices k in fip gaugino/gluino decays refer to the ordering of the indices 
in the fip couplings in the superpotential. The convention is as in ref. [42]. 

Next a number of parameters derived from the SUSY lagrangian must be given. 
These are: the ratio of Higgs VEVs, tan/?, and the scalar Higgs mixing angle, a; the 
mixing parameters for the Higgses, gauginos and the sleptons; the trilinear couplings; 
and the Higgsino mass parameter fx. 

Finally the logical variable RPARTY must be set . FALSE . if R-parity is violated, 
and the fip couplings must also be supplied; otherwise not. 

Details of the FORMAT statements employed can be found by examining the sub- 
routine HWISSP. 

HWISSP reads the data from UNIT = LRSUSY (defauh LRSUSY = 66). If the data 
are stored in a fort. LRSUSY file on a UNIX system^ no further action is required, 
but if the data are to be read from a file named f name . dat then appropriate OPEN 
and CLOSE statements must be added by hand: 

OPEN (UNIT=LRSUSY , FDRM= ' FORMATTED ' , STATUS= ' UNKNOWN ' , FILE= ' f name . dat ' ) 
CALL HWISSP CLOSE (UNIT=LRSUSY) 

A number of sets of SUSY parameter files, produced using ISAJET, for the stan- 
dard LHC SUGRA and GMSB points are available from the ISAWIG home page: 
http : / /www . hep . phy . cam . ac . ukTh ^richardn/HERWIG/ISAWIG/ 

^Or the equivalent, e.g. fortLRSUSY.dat on a VAX system. 
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Decaying 
Particle 


No. of 
products 


Type of mode 


Order of decay products 


1st 


2nd 


3rd 


top 


2 


two-body to Higgs 


Higgs 


bottom 




3 


three-body via Higgs/W 


quarks or leptons 
from W/Higgs 


bottom 


gluinos 


2 


without gluon 


any order 




with gluon 


gluon 


colour 
neutral 




3 


R-parity conserved 


colour 
neutral 


quark or antiquark 


squark 
or slepton 


2 


gaugino/gluino 
quark/lepton 


gaugino 
gluino 


quark 
lepton 




3 


weak 


sparticle 


particles from 
W decay 


squarks 


2 


lepton number violated 


quark 


lepton 




baryon number violated 


quark 


quark 




sleptons 


2 


lepton number violated 


quark or antiquark 




Higgs 


2 


(s)quark-anti(s)quark 


(s)quark or anti(s)quark 




(s) lepton- anti (s) lepton 


(s)lepton or anti(s)lepton 




3 


all three-body 


colour 
neutral 


quark or antiquark 
lepton or antilepton 


gaugino 


2 


squark-quark 


quar. 


-c or squark 




slepton-lepton 


lepton or slepton 




3 


R-parity conserved 


colour 
neutral 


quark or antiquark 
lepton or antilepton 


gaugino 
or gluino 


3 


R-parity violating 


particles in the order k as 
in the superpotential 



Table 4: SUSY decay product ordering. 



3.5.2 Processes 



The implementation of supersymmetric particle production processes in lepton-anti- 
lepton and hadron-hadron collisions is described in sections |4.4| and [4.7| , respectively. 
As in SM processes, we do not include any higher-order QCD corrections to the 
relevant matrix elements, even though these are now known in many cases. This is 
in order to avoid double counting of corrections generated approximately by parton 
showering. 

The procedure by which the MSSM matrix elements describing the elementary 
hard subprocesses are interfaced to the initial- and final-state parton showers is 



similar to that described in sections and |3.2| for the SM case. However, at present 
showers are generated from partons but not from spartons, whose short lifetimes 
should make this a reasonable approximation. 
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We include matrix elements and spin correlations in all SUSY decays as described 
in [21]. Breit-Wigner mass distributions have been implemented for all unstable 
SUSY particles, according to their widths as given in the input file. For a list of the 
most relevant decay modes of MSSM particles, see ref. [22]. 

One difference between the SM and MSSM implementations of the Higgs decay 
channels should be mentioned. Whereas for the SM Higgs boson all decay rates are 
calculated by HERWIG itself, in terms of the other (known) parameters of the SM, 
in the case of the MSSM scalars these are passed to the generator through the data 
files. 

A more detailed description of the new SUSY reactions in HERWIG, along with 
the relevant formulae for the hard scattering processes involved, can be found in [22]. 

We include the possibility of violating R-parity both in the decays of the sparti- 
cles and in the initial hard subprocess. The flp model we consider is specified by a 
spectrum which can be given in either the SUGRA, GMSB or a more general MSSM 
scenario, and a set of couplings at the weak scale. We include only the tri-linear 
couplings and not the bi-linear terms which mix the leptons and gauginos; a recent 
review of these models can be found in [43]. However we do include the possibility 
that more than one of the fip couplings can be non-zero. All the two-body squark and 
slepton and three-body gaugino and gluino decay modes, and resonant production 
processes in hadron-hadron collisions are included, as well as a range of production 
processes in e'^e~ collisions. A decay matrix element is implemented for the three- 
body decays. The colour structure of these events is very different from that of the 
MSSM [42], due to the presence of the baryon- number violating (p) vertex. This 
means that a new subroutine HWBRCN is required to handle the colour connections be- 
tween jets in this case. A full discussion of the colour connections in these processes 
and the matrix elements for the cross sections and decays can be found in [42]. 

3.5.3 Related changes 

A large number of changes have been made in HERWIG to enable SUSY processes 
to be included in hadron-hadron collisions. The main changes are: 

• The subroutine HWDHQK has been replaced by HWDHOB which does both heavy 
quark and SUSY particle decays. 

• The subroutines HWBCON, HWCGSP and HWCFOR have been adapted to handle the 
colour connections found in normal SUSY decays. 

• The subroutine HWBRCN has been included to deal with the inter-jet colour 
connections arising in fip SUSY. Also HWCBVI, HWCBVT and HWCBCT have been 
added to handle the hadronization of baryon number violating (p) SUSY decays 
and processes. If the variable RPARTY= . TRUE . (default) then the old HWBCON 
colour connection code is used, else the new HWBRCN. 
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3.6 Spin correlations 

In addition to the spin correlations in parton showers and vector boson decays dis- 
cussed above, in versions 6.4 and higher spin correlation effects are included in pro- 
cesses where top quarlcs, r leptons and SUSY particles are produced, as described 
in [21]. At the moment the effects are only calculated for the production of these 
particles in the following processes: IPR0C=100-199, 700-799, 800-899, 1300-1399, 
1400-1499, 1500-1599, 1700-1799, 2000-2099, 2800-2825, 3000-3030, 4000-4199. How- 
ever, if these particles are produced in other processes, the spin correlation algorithm 
will still be used to perform their decays. The correlations are also included for the 
decay of the MSSM Higgs bosons, regardless of how they are produced. 

The spin correlations are controlled by the logical variable^ SYSPIN [.TRUE.] 
which switches the correlations on. If required the correlations are initialised by the 
new routine HWISPN. This routine initialises the two, three and four body matrix 
elements. 

The three and four body matrix elements can be used separately to generate the 
decay distributions without spin correlation effects. These are switched on by the 
switches THREEB [ . TRUE . ] for three body decay and FOURB [ . FALSE . ] for four body 
decays. The four body decays are only important in SUSY Higgs studies, and have 
small branching ratios. However, they take some time to initialise and are therefore 
switched off by default. 

The initialisation of the spin correlations and/or decay matrix elements can be 
time consuming and we have therefore included an option to read/write the informa- 
tion. The information is written to unit LWDEC [88] and read from LRDEC [0] . If either 
are zero the data is not written/read. If IPRINT=:2 then information on the branch- 
ing ratios for the decay modes and the maximum weights for the matrix elements is 
outputted. 

If the spin correlation (SYSPIN) or matrix element switches (THREEB, FOURB) are 
. TRUE . , then the matrix element codes (NME entries) for the decays concerned are 
not used; the calculated matrix elements are used instead. 

When we included spin correlations in HERWIG6.4 [21] we did not include either 
R-parity violating decays or decays producing gravitinos in the algorithm. This 
led to HERWIG stopping when such decays were included. This of course could be 
stopped by switching the spin correlations off, i.e. SYSPIN= . FALSE . . In version 6.5, 
we have included the relevant matrix elements for R-parity violating decays and hard 
processes and decays producing gravitinos. At the same time we have made changes 
so that at both the initialisation and event generation stages many of the terminal 
warnings which were caused by the code not having the correct matrix elements are 
now information-only warnings. If you still get terminal error messages from any of 
the spin correlation routines please let us know. 

^Default values for input variables are shown in square brackets. 
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The effect of polarization for incoming leptonic beams in MSSM and SUSY 
processes lias also been included. These effects are included both in the production 
of SUSY particles and via the spin correlation algorithm in their decays. 

3.7 Hadronization 

For a general hard process in hadron-hadron collisions, we have to consider: (a) the 
representation of the incoming partons as constituents of the incident hadrons; (6) the 
conversion of the emitted partons into outgoing hadrons; (c) the 'underlying soft 
event' associated with the presence of spectator partons. 

The first of these is dealt with through the use of non-perturbative parton dis- 
tribution functions, which are discussed below in section [4.1.1| , and by the remnant 
hadronization model. The cluster model for hadron formation, remnant hadroniza- 
tion and the underlying event is as follows. 

3.7.1 Cluster model 

The preconfinement property mentioned in section |1| is used by HERWIG as the basis 
for a simple hadronization model which is local in colour and independent of the 
hard process and the energy [4,7]. 

After the perturbative parton showering, all outgoing gluons are split non-pertur- 
batively, into light quark- ant iquark or diquark-antidiquark pairs (the default option 
is to disallow diquark splitting). At this point, each jet consists of a set of outgoing 
quarks and antiquarks (also possibly some diquarks and antidiquarks) and, in the 
case of spacelike jets, a single incoming valence quark or antiquark. The latter 
is replaced by an outgoing spectator carrying the opposite colour and the residual 
flavour and momentum of the corresponding beam hadron. 

In the limit of a large number of colours, each final-state colour line can now 
be followed from a quark/anti-diquark to an antiquark/diquark with which it can 
form a colour- singlet cluster.^ By virtue of pre-confinement, these clusters have a 
distribution of mass and spatial size that peaks at low values, falls rapidly for large 
cluster masses and sizes, and is asymptotically independent of the hard subprocess 
type and scale. 

The clusters thus formed are fragmented into hadrons. If a cluster is too light 
to decay into two hadrons, it is taken to represent the lightest single hadron of its 
flavour. Its mass is shifted to the appropriate value by an exchange of 4-momentum 
with a neighbouring cluster in the jet. Similarly, any diquark-antidiquark clusters 
with masses below threshold for decay into a baryon-antibaryon pair are shifted to 
the threshold via a transfer of 4-momentum to a neighbouring cluster. 

^The situation when baryon number is violated is more complicated and is discussed in [44] for 
the Standard Model and in [42] for R-parity violating SUSY models. 
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Those clusters massive enough to decay into two hadrons, but below a fission 
threshold to be specified below, decay isotropically into pairs of hadrons selected 
in the following way. A fiavour / is chosen at random from among u, d, s, the six 
corresponding diquark fiavour combinations, and c. For a cluster of fiavour /1/2, 
this specifies the fiavours fif and //2 of the decay products, which are then selected 
at random from tables of hadrons of those fiavours. See section |^ for details of the 
hadrons included. The selected choice of decay products is accepted in proportion to 
the density of states (phase space times spin degeneracy) for that channel. Otherwise, 
/ is rejected and the procedure is repeated. 

The above method of selection for cluster decays is simple and fast but does not 
automatically satisfy constraints such as strong isospin symmetry. The decay rate 
into hadrons of a certain fiavour depends on the average phase space for channels 
involving that fiavour. Thus, for example, the existence of the 77 or 77', with the same 
quark content as the 7r°, leads to a slight reduction of direct vr" production relative to 
TT~^ and TT^. Quantitatively, the effect is too small to be observed even with the high 
statistics of the LEPl data. However, the method can give rise to strange effects 
if the particle data tables are extended, and modifications to avoid this have been 
proposed [45]. 

In the decays of clusters to 77 or 77', the parameter ETAMIX gives the fjs/vo mixing 
angle in degrees (default = -20). Rates are not very sensitive to its exact value, as 
the 7]'/ 1] suppression is dominated by mass effects in the cluster model. See section |^ 
for more details. 

A fraction of clusters have masses too high for isotropic two-body decay to be a 
reasonable ansatz, even though the cluster mass spectrum falls rapidly (faster than 
any power) at high masses. These are fragmented using an iterative fission model 
until the masses of the fission products fall below the fission threshold. In the fission 
model the produced fiavour / is limited to u, d or s and the product clusters fif and 
7/2 move in the directions of the original constituents /i and /2 in their cm. frame. 
Thus the fission mechanism is not unlike string fragmentation [46]. 

In HERWIG there are three main fission parameters, CLMAX, CLPOW and PSPLT. 
The maximum cluster mass parameter CLMAX and CLPOW specify the fission threshold 
Mf according to the formula 

^CLPDw ^ CLMAX^^P™ + (mi + 7712)'^'™ , 

where tt^i and m2 are the quark mass parameters RMASS(i) for fiavours /i and /2 
(see section |3.2| ). The parameter PSPLT specifies the mass spectrum of the produced 
clusters, which is taken to be M^^^^"^ within the allowed phase space. Provided the 
parameter CLMAX is not chosen too small (the default value is 3.35 GeV), the gross 
features of events are insensitive to the details of the fission model, since only a small 

"^Except for those containing a 'perturbative' quark when CLDIR =1 — see below. 
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fraction of clusters undergo fission. However, the production rates of high-p^ or heavy 
particles (especially baryons) are affected, because they are sensitive to the tail of 
the cluster mass distribution. The default value of the power CLPOW is 2. Smaller 
values will increase the yield of heavier clusters (and hence of baryons) for heavy 
quarks, without affecting light quarks much. For example, the default value gives no 
6-baryons (for the default value of CLMAX) whereas CLPOW =1.0 makes the ratio of 
6-baryons to 6-hadrons about 1/4. 

There is also a switch CLDIR for cluster decays. If CLDIR = 1 (the default) then 
a cluster that contains a 'perturbative' quark, i.e. one coming from the perturbative 
stage of the event (the hard process or perturbative gluon splitting) 'remembers' its 
direction. Thus when the cluster decays, the hadron carrying its flavour continues in 
the same direction (in the cluster cm. frame) as the quark. This considerably hard- 
ens the spectrum of heavy hadrons, particularly of c- and 6-flavoured hadrons. It also 
introduces a tendency for baryon-antibaryon pairs preferentially to align themselves 
with the event axis (the 'TPC/27 string effect' [47]). CLDIR = turns off this op- 
tion, treating clusters containing quarks of perturbative and non-perturbative origin 
equivalently. In the CLDIR = 1 option, the parameter CLSMR (default = 0.0) allows 
for a gaussian smearing of the direction of the perturbative quark's momentum. The 
smearing is actually exponential in 1 — cos 9 with mean CLSMR. Thus increasing CLSMR 
decorrelates the cluster decay from the initial quark direction. 

The process of 6-quark hadronization requires special treatment and the results 
obtained using HERWIG are still not fully satisfactory. Generally speaking, it is 
difficult to obtain a sufficiently hard B-hadron spectrum and the observed 6-meson/6- 
baryon ratio. These depend not only on the perturbative subprocess and parton 
shower but also on non-perturbative issues such as the fraction of 6-flavoured clusters 
that become a single B meson, the fractions that decay into a B meson and another 
meson, or into a 6-baryon and an antibaryon, and the fraction that are split into more 
clusters. Thus the properties of 6-jets depend on the parameters RMASS(5), CLMAX, 
CLPOW and PSPLT in a rather comphcated way. In practice these parameters are tuned 
to global final-state properties and one needs extra parameters to describe 6-jets. 

A parameter BILIM has therefore been introduced to allow clusters somewhat 
above the Btt threshold mass Mth to form a single B meson if 

M < Mii^ = (1 + BlLIM)Mth . 

The probability of such single-meson clustering is assumed to decrease linearly for 
Mth < M < Miim. This has the effect of hardening the B spectrum if BILIM is 
increased from the default value of zero. In addition, in version 6, the parameters 
PSPLT, CLDIR and CLSMR have been converted into two-dimensional arrays, with the 
first element controlling clusters that do not contain a fe-quark and the second those 
that do. Thus tuning of 6-fragmentation can now be performed separately from other 
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flavours, by setting CLDIR(2) = 1 and varying PSPLT(2) and CLSMR(2). By reducing 
the value of PSPLT(2), further hardening of the B-hadron spectrum can be achieved. 

3.7.2 Underlying soft event 

In hadron-hadron and lepton-hadron collisions there are 'beam clusters' containing 
the spectators from the incoming hadrons. In the formation of beam clusters, the 
colour connection between the spectators and the initial-state parton showers is cut 
by the forced emission of a soft quark-antiquark pair. The underlying soft event in a 
hard hadron-hadron collision is then assumed to be a soft collision between these two 
beam clusters. In a lepton-hadron collision the corresponding 'soft hadronic remnant' 
is represented by a soft collision between the beam cluster and the adjacent cluster, 
i.e. the one produced by the forced emission mentioned above. 

The model used for the underlying event is based on the minimum-bias pp event 
generator of the UA5 Collaboration [48], modified to make use of our cluster frag- 
mentation algorithm. This model is explained in the following subsection. 

Adding 10000 to the HERWIG process code IPROC suppresses the underlying 
event, in which case the beam clusters are simply fragmented like other clusters, 
without any soft collision. The parameter PRSOF enables one to produce an underly- 
ing event in only a fraction PRSOF of events (default = 1.0). Adding 10000 to IPROC 
is thus equivalent to setting PRSOF = 0. 

A parameter BTCLM is available to users to adjust the mass parameter equivalent 
to PMBMl (see below) in remnant cluster formation. Its default value, 1.0, is identical 
to previous versions. There is also an option for the special treatment of the splitting 
of clusters containing hadron (or photon) remnants. lOPREM = gives the fragments 
a gaussian mass spectrum typical of soft processes. When lOPREM = 1 (default), the 
child containing the remnant is treated as before but the other cluster, containing a 
perturbative parton, is treated as a normal cluster, with mass spectrum M^^^^^ . 

Two special remnant 'particles' have been defined: 'REMG ' with IDHW = 71, 
IDHEP = 98 and 'REMN ' with IDHW = 72, IDHEP = 99. These are remnant pho- 
tons and nucleons respectively. They are identical to photons and nucleons, except 
that gluons are labelled as valence partons and, for the nucleon, valence quark dis- 
tributions are set to zero. They are used by an external package for simulating 
multi-parton interactions, called JIMMY [49]. See section for further details. 

3.7.3 Minimum bias processes 

The minimum-bias event generator of the UA5 Collaboration [48] starts from a 
parametrization of the pp inelastic charged multiplicity distribution as a negative 
binomial distribution. In HERWIG version 6, the relevant parameters are made avail- 
able to the user for modification, the default values being the UA5 ones as used in 
previous versions. These parameters are given in table |^. 
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Name 


Description 


Default 


PMBNl 


a in n = as^ + c 


9.110 


PMBN2 


b in n — as'' + c 


0.115 




in /^ 1 


Q "^nn 
— y .ouu 


PMBKl 


0, in 1 / /c = a In s + 6 


0.029 


PMRK9 


n m X 1 ic — n In c —1— n 
\j 111 J- 1 th — (X ill o 


—0 104 


PMBMl 


a in (M — mi — m2 — a)e~^^ 


0.4 


PMBM2 


6 in (M — mi — m2 — a)e^^*^ 


2.0 


PMBPl 


Pt slope for w 


5.2 


PMBP2 


Pt slope for s, c 


3.0 


PMBP3 


Pt slope for qq 


5.2 



Table 5: Soft/min.bias parameters. 



The first three parameters control the mean charged multiplicity n at cm. energy 
as indicated. The next two specify the parameter k in the negative binomial 
charged multiplicity distribution, 

^ ^ n\T{k) (l+n/fc)"+*^' 

The parameters PMBMl and PMBM2 describe the distribution of cluster masses M in 
the soft collision. These soft clusters are generated using a fiat rapidity distribution 
with gaussian shoulders. The transverse momentum distribution of soft clusters has 
the form 

P(pt) ocp^exp + M^) , 

where the slope parameter b depends as indicated on the flavour of the quark or 
diquark pair created when the cluster was produced. As an option, for underlying 
events the value of ^/s used to choose the multiplicity n may be increased by a 
factor ENSOF to allow for an enhanced underlying activity in hard events. The actual 
charged multiplicity is taken to be n plus the sum of the moduli of the charges of 
the colliding hadrons or clusters. 

There is now also an interface to the multiple-interaction model JIMMY [49]. 
For this purpose, several routines have been added or modified. New are HWHREM for 
identifying and cleaning up the beam remnants and HWHSCT to administer the extra 
scatters. Minor modifications to HWBGEN and HWSBRN suppress energy conservation 
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errors when ISLENT = — 1; HWSSPC has an improved approximation for remnant mass 
at high energies; and HWUPCM improves safety against negative square roots. 



3.8 Spacetime structure 

The space-time structure of events is available for all types of subprocess. The pro- 
duction vertex of each parton, cluster, unstable resonance and final-state particle is 
supplied in the VHEP array of /HEPEVT/. Set PRVTX = .TRUE, to include this infor- 
mation when printing the event record (120 column format). The units are: x,y,z 
in mm and t in mm/c. In the case of partons and clusters the production points are 
always given in a local coordinate system with its origin at the relevant hard subpro- 
cess. This helps to separate the fermi-scale partonic showers from millimetre-scale 
distances possible in particle decays, for example the partonic decays of heavy (c, b) 
hadrons. The vertices of hadrons produced in cluster decays are always corrected 
back into the laboratory coordinate system. 

It is possible to vary the principal interaction point, assigned to the cm. frame 
entry in /HEPEVT/ (with ISTHEP = 103), by setting PIPSMR = .TRUE. The smearing is 
generated by the routine HWRPIP according to a triple gaussian given by parameters 
VIPWID(I) (I = 1,2,3 for x,y,z widths): the default values correspond to LEPl. 

It is also possible to veto particle decays that would occur outside a specified 
volume by setting MAXDKL = .TRUE. Each putative decay is tested in HWDXLM and if 
the particle would have decayed outside the chosen volume it is frozen and labelled 
as final state. Using lOPDKL = 1,2 selects a cylindrical or spherical allowed region 
(centred about the origin): then parameters DXRCYL, DXZMAX or DXRSPH specify the 
dimensions of the region. 

3.8.1 Particle decays 

Lepton and hadron lifetimes (in seconds) are supplied in the array RLTIM. In the 
case of MSSM (s)particles, including Higgs states, RLTIM values are entered through 
the input files (see discussion in section |3.5.2 ). The lifetimes of heavy quarks (top 



and any hypothetical extra generations) and weak bosons (including the SM Higgs) 
are derived from their calculated or specified widths in HWUDKS, whilst light quarks 
and gluons are given an effective minimum width that acts as a lifetime cutoff — see 
below. All particles whose lifetimes are larger than PLTCUT are set stable. 

The proper (i.e. rest-frame) time t* at which an unstable lepton or hadron decays 
is generated according to the exponential decay law with mean lifetime (t*) = r = 
RLTIM: 

Prob(proper time > t*) = exp(— t*/r) . 
The laboratory-frame decay time t and distance travelled d are obtained by applying 



a boost: t = 'yt*, d = f^'yt* where (3 = v/c and 7 = l/yl — (3"^. The production 
vertices of the daughter particles are then calculated by adding the distance travelled 
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by the mother particle as given above to its production vertex. The mean hfetime r 
of a particle is set, taking into account its width and virtuality, by: 



This formula is used for all particles: light partons; heavy quarks and weak bosons, 
which have appreciable widths; resonances; and unstable leptons. It interpolates 
between r = h/V for a particle that is on mass-shell and r = fi\f(fl[(f' — M^) for 
one that is far off mass-shell. 

3.8.2 Parton showers 

The above prescription, based on an exponential proper lifetime distribution, is also 
used to describe parton showers. For light quarks and gluons, whose natural widths 
are small, this could lead to unreasonably large distances being generated in the final, 
low virtuality steps of showering. To avoid this they are given a width F = VMIN2/M; 
the parameter VMIN2 (default value 0.1 GeV^) acts effectively as a lower limit on a 
parton's virtuality. This is particularly important for the forced splitting of gluons 
(see section |]7[l|), which uses r = ;iRMASS(l3)/VMIN2. 

3.8.3 Hadronization 

In the case of a cluster its initial production vertex is taken as the midpoint of a 
line perpendicular to the cluster's direction of travel and with its two ends on the 
trajectories of the constituent quark-antiquark pair. If such a cluster undergoes a 
forced sphtting to two clusters the string picture is adopted. The vertex of the light 
quark pair is positioned so that the masses of the two daughter clusters would be the 
same as those for two equivalent string fragments. The production vertices of the 
daughter clusters are given by the first crossing of their constituent qq pairs. The 
production positions of primary hadrons from cluster decays are smeared, relative to 
the cluster position, according to a gaussian distribution of width 1/ (cluster mass). 

3.8.4 Colour rearrangement 

HERWIG version 6 contains a colour rearrangement model based on the space-time 
structure of an event at the end of the parton shower. This is illustrated in the simple 
example shown below where showering results in a colour-neutral qggq final state. In 
the conventional HERWIG hadronization model (corresponding to the default value 
of the reconnection parameter, CLRECO = .FALSE.), after a non-perturbative splitting 
of the final-state gluons, colour singlet clusters are formed from neighbouring qq 
pairs: {ij){pq){kl). However when CLRECO = .TRUE, the program first creates colour 
singlet clusters as normal but then checks all (non-neighbouring) pairs of clusters 
to test if a colour rearrangement lowers the sum of the clusters' spatial sizes added 




v/(g2 - M2)2 + (FgVM)2 ■ 
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in quadrature. A cluster's size dij is defined to be the Lor entz- invariant space-time 
distance between the production points of its constituent quark qi and antiquark qj. 
If an allowed alternative is found, that is, {ij){kl) — > {il){jk) such that + > 
Mid^ + Mfe/P) then it is accepted with a probability given by the parameter PRECO 
(default value 1/9). 




Note that not all colour rearrangements are allowed, for instance in the example 
shown {ij){pq) {i(l){jp) is forbidden since the cluster {jp) is a colour octet — it 
contains both products from a non-perturbative gluon splitting. 

Multiple colour rearrangements are considered by the program, as are those be- 
tween clusters in jets arising from a single, colour neutral source, for example 
decay (as shown above), or due to more than one source, for example e+e" — > 
W'^W~ 4 jets. In the latter case a new parameter, EXAG, is available to exag- 
gerate the lifetime of the or any other weak boson, so that any dependence of 
rearrangement effects on source separation can be investigated. 

The CLRECO option can be used for all the processes available in HERWIG. Note, 
however, that before using the program with CLRECO = .TRUE, for detailed physics 
analyses the default parameters should be retuned to LEP data with this option 
switched on. 



3.8.5 B-B mixing 

When MIXING = .TRUE., particle- antiparticle mixing for mesons is implemented. 
The probability that a meson is mixed when it decays is given in terms of its lab-frame 
decay time t by: 

1 cos{Xtm/cTE) 
Pmix(i) = 2 ~ 2cosh{Ytm/crE) ' 

where X = AM/T, Y = Ar/2r and m, r, E are the B° mass, lifetime and energy. 
The ratios X and Y are stored in XMIX(I) and YMIX(I), I = 1,2 for g = s,d. 
Whenever a neutral B meson occurs in an event, a copy of the original entry is 
always added to the event record, with ISTHEP = 200, which gives the particle's 
flavour at the production (or cluster decay) time. This is in addition to the usual 
decaying particle entry with ISTHEP = 199. 
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Name 


Description 


Default 


PARTI 


Type of particle in beam 1 


'PBAR ' 


PART2 


Type of particle in beam 2 


'P ' 


PBEAMl 


Momentum of beam 1 


900.0 


PBEAM2 


Momentum of beam 2 


900.0 


IPROC 


Type of process to generate 


1500 


MAXEV 


Number of events to generate 


100 



Table 6: Main program variables. 



4. Processes 
4.1 Beams 

As indicated in table ^, a number of variables must be set in the main program 
HWIGPR to specify what is to be simulated. The beam particle types PARTI, PART2 
can take any of the values NAME listed in table |^. 
In the case of point-like photon/QCD pro- 
cesses, IPROC = 5000-5999, the first particle 
must be the photon or a lepton. In addition, 
beams 'K+ ' and 'K- ' are supported for 
minimum bias non-diffractive soft hadronic 
events (IPROC = 8000) only. In the case that 
the beam momenta PBEAMl and PBEAM2 are 
not equal, the default procedure (USECMF = 
.TRUE.) is to generate events in the beam- 
beam centre-of-mass frame and boost them 
back to the laboratory frame afterwards. 

In hadronic processes with lepton beams 
(e.g. photoproduction in ep), the lepton —> 
lepton + photon vertex uses the full transverse-momentum dependent splitting func- 
tion, with exact light-cone kinematics, i.e the Equivalent Photon Approximation 
(EPA). This means that the photon-hadron collision has a transverse momentum in 
the lepton-hadron frame and must be boosted to a frame where it has no transverse 
momentum. Thus the c.m.f. boost described above is always used in these processes, 
regardless of the value of USECMF. The correct lower energy cutoff appropriate to the 
hadronic process is applied to the photon. The of the photon is generated within 
the kinematically allowed limits, or the user-defined limits Q2WWMN and Q2WWMX (de- 
faults and 4) whichever is more restrictive.* Similarly for the photon's light-cone 
momentum fraction, with user-defined limits YWWMIN and YWWMAX (default and 1). 

*The WW in parameter names is a relic from earlier versions that used the less accurate Weizsacker- 
Williams approximation. 





NAME 




NAME 


e+ 


'E+ ' 


e 


'E- ' 




'MU+ ' 




'MU- ' 




'NU_E ' 




'NU_EBAR ' 




'NU_MU ' 




'NU_MUBAR' 




'NU_TAU ' 




'NU_TAUBR' 


V 




P 


'PBAR ' 


n 




n 


'NBAR ' 


7r+ 


'PI+ ' 




'PI- ' 


7 


'GAMMA ' 







Table 7: Beam particles. 
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Together with the Bjorken y-variable hmits YBMIN and YBMAX, this allows different 
ranges for the tagged and untagged photons in two-photon DIS. 

4.1.1 Parton distributions 

The parton momentum fraction distributions of the beam particles are used in the 
generation of initial-state parton showers and also in the non-perturbative process 
of linking the shower with the beam hadron and its remnant. Since the parton 
showering is done in leading-logarithmic order, there is no strong motivation to use 
next-to-leading order parton distributions, although this has become customary since 
the most up-to-date distributions are deduced from next-to-leading order fits to (in- 
clusive) data. Thus the most common option is to use the interface to the PDFLIB 
parton distribution library [50]. 

The HERWIG interface is compatible with PDFLIB version 4. AUTPDF should be 
set to the author group as listed in the PDFLIB manual, e.g. 'MRS' , and MODPDF to the 
set number in the new convention. It is permissible to choose the PDFLIB set inde- 
pendently for each of the two beams. For example, to use MRS D- for the proton and 
Gordon- Storrow set 1 for the photon in 7-hadron or lepton-hadron collisions, one sets: 

AUTPDF (2)= 'MRS' 
MODPDF (2) =28 
AUTPDF(1)='GS' 
MODPDF (1) =2 

If the PDFLIB interface is not used, the parton distributions are chosen from the 
HERWIG internal sets according to the value of the parameter NSTRU. The default 
parton distributions in HERWIG versions prior to 6.3 were very old and did not in- 
clude fits to any of the HERA data. Therefore several new PDFs have been included 
in versions 6.3 and higher. These are shown in table 

It should be noted that we have only added leading-order fits because the evo- 
lution algorithms in HERWIG, in particular the backward-evolution algorithm for 
initial-state parton showering, are only leading-order and therefore inconsistencies 
could occur with next-to-leading-order distributions. 

The new default structure function set NSTRU=8 is the average of two of the 
published fits [51], because this has been found [52] to be closer to the central value 
of more recent next-to-leading-order fits. The other fits can then be used to assess 
the effects of varying the high-x gluon. 

These new HERWIG parton distributions are only available for nucleons. For 
pion beams, either the old NSTRU=1,2 pion sets or PDFLIB should be used. 

For photons, the default is to use the Drees-Grassie parton distributions [53]. 
The heavy quark content of the photon uses the corrections to the Drees-Grassie 
distribution functions for light quarks, calculated by Drees and Kim [54]. There is 
also an interface to the Schuler-Sjostrand [55] parton distribution functions for the 
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NSTRU 


Description 


6 
7 

8 


Central as and gluon leading-order fit of [51] 
Higher gluon leading-order fit of [51] 

Average of central and higher gluon leading-order fits of [51] 



Table 8: New internal MRST parton distributions. 



photon, version 2. These appear as PDFLIB sets with author group 'SaSph', but are 
actually implemented via a call to their SASGAM code. The value in MODPDF specifies 
the set (1-4 for ID [recommended set],lM, 2D,2M), whether the Bethe-Heitler process 
is used for heavy flavours (add 10), whether the P^-dependence is included (add 20), 
and which of their models is used (add 100 times their IP2 parameter). 

An option to damp the parton distributions of off mass-shell photons relative 
to on-shell photons, according to the scheme of Drees and Godbole [56] has been 
introduced. The adjustable parameter PHOMAS defines the crossover from the non- 
suppressed to suppressed regimes. Recommended values lie in the range from QCDLAM 
to 1 GeV. The default value PHOMAS = corresponds to no suppression, as in previous 
versions. 



4.2 Summary of subprocesses 



We give in table ^ a list of the currently available hard subprocesses IPROC. More 
detailed descriptions are given in sections [4.3H4.12| , and then in section |4.13| there 
are instructions to users on how to add a new process. 



IPROC 



Process 



100 
100+IQ 

107 

110 
110+IQ 

120 
120+IQ 

127 
150+IL 



qq{g) (all q flavours) 
£+£- qq{g) (IQ = 1, 2, 3, 4, 5, 6 for g = d, m, s, c, 6, t) 
i'^i^ gg{g) (fictitious process) 

qqg (all flavours) 

qqg (IQ as above) 
i~^£^ —>■ qq (all flavours, no hard gluon correction) 

gg (IQ as above, no hard gluon correction) 

gg (fictitious process, no hard gluon correction) 
£+£- i'l' (IL =1,2,3 for f = e, fi, r, N.B. i ^ i') 



200 
250 



i+t 



W~^W (see sect. [4.3.2| on control of W/Z decays) 
(see sect. [4.3.2| on control of W/Z decays) 



300 
300+IQ 
306+IL 
310, 311 



£+t 
£+t 



70770 
7O rrO 
7O rrO 



Z^qq (all flavours) 
Z^qq (IQ as above) 
Z^U (IL as above) 

z^w+w-, Z^Z^Z^ 

Table ^: Process codes. (Continues) 
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IPROC 


Process 


312 
399 


i+t 
£+t 


ZOiJgV ^ anything 


400+ID 


£+£- 


i^uH^M + ^^^~HsM (ID as in IPROC = 300 + ID) 


500+ID 
550+ID 


(ID= 

£+£- 


£+£"77 ^ £+£-qq/££/W-^W- 
=0-10 as in IPROC = 300 + ID) 

£ue'yW £vm' / ££' (lD=0-9 as in IPROC = 300 + ID) 


600 
600+IQ 


qqgg, qqq'q' (all q flavours) 
qqgg, qqq'q' (IQ as above) 
After generation, IHPRO is subprocess (see sect. [4.3.5|) 


700-99 
700 
710 

706+4IN1+IN2 
730 

728+2IC1+IC2 

740 
736+5IL 
737+5IL 
738+5IL 
739-I-5IL 
740+5IL 

760 
757+4IQ 
758+4IQ 
759+4IQ 
760+4IQ 


Minimal Supersymmetric Standard Model (MSSM) processes 

2-sparticle processes (sum of 710, 730, 740 and 760) 
t^£~ —>■ neutralino pairs (all neutralinos) 

XiNiXiN2 (INI ,2=neutralino mass eigenstate) 
—>■ chargino pairs (all charginos) 

XiciXic2 (ICl , 2=chargino mass eigenstate) 
£~^£~ —>■ slepton pairs (all flavours) 
i^-r ^ IlJU (IL = 1,2,3 for i = e,fi, f) 
i^r ejl (IL as above) 
£+£- IlIr (IL as above) 

Ir£*r (IL as above) 
£+£- vlvI (IL = 1,2,3 for v^, u^, u^) 
£'^£^ squark pairs (all flavours) 

^ qL,RqlR (IQ = 1, 2, 3, 4, 5, 6 for g = d, u, s, c, 6, i) 
— > qiql (IQ as above) 
£~^£~ — > qiqii (IQ as above) 
£~^£~ — > qaqji (IQ as above) 


800-99 

800 

810 
810+IN 

820 
820+IC 

830 

840 

850 

860 

870 

867-F3IL1+IL2 


R-parity violating supersymmetric processes 

Single sparticle production, sum of 810-840 

£~^£~ — >■ x^Ui, (all neutralinos) 

£~^£~ — > y?M^',-, f IN=neutralino mass state) 

£'^£^ X'ej' (all charginos) 

£'^£^ Xic^ti (lC=chargino mass state) 

£+£- UiZ^ and IfW- 

£+£- UihynyA^ and £+£- IfH- 

£+£' ua 

Sum of 870 and 880 

£+£- £+£-, via LIE only 

£+£- Ali^ils (ILI, 2=1,2,3 for e,/i,r) 

Table ^: Process codes. (Continues) 
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IPROC 

880 

^77+3lQl+IQ2 



Process 



dd, via LLE and LQD 

c?iLi^iL2 (IQI, 2=1,2,3 for d,s,b) 



910 
920 
960 
970 
955 
965 
965 



£+t 



H+H 



1000+ID 



ttH^^ (ID as in IPR0C=300+ID) 



1110+IQ 
1116+IL 
1120+IQ 
1126+IL 
1130+IQ 
1136+IL 

1140 

1141 

1142 

1143 

1144 

1145 



i+t 
i+t 
i+t 
i+t 



qqh^ (IQ as in IPR0C=100+IQ) 

(IL=1,2,3 for e,/i,r) 
qqH^ (IQ as in IPR0C=100+IQ) 

(IL=1,2,3 for e,/i,r) 
qqA^ (IQ as in IPR0C=100+IQ) 
£+£-^0 (IL=1,2,3 for e,/i,r) 
duH^+ ch. conj. 
scH~^+ ch. conj. 
biH^+ ch. conj. 
e UeH~^+ ch. conj. 
/i u^H^+ ch. conj. 
tUt-H^+ ch. conj. 



1200-99 



Reserved for other processes 



1300 
1300+IQ 

1350 
1350+IL 

1399 



qq 
qq 
qq 
qq 
qq 



zVi 



q'q' (all flavours) 
q'q' (IQ = 1,2,3,4,5,6 for q 
U (all lepton species) 
U (IL = 1 — 6 for £ = e, i^e, /W, z/^ 
anything 



(i, M, s, c, 6, t) 



^l1 



etc. 



1400 
1400+IQ 

1450 
1450+IL 

1499 



qq 
qq 
qq 
qq 
qq 



q'q" (all flavours) 

q'q" {q' or q" as above) 

iug (all lepton species) 

£ue (IL = 1,2,3 for i = e,fi,T) 

anything 



1500 



QCD 2 — > 2 hard parton scattering 

After generation, IHPRO is subprocess (see sect. [4.6.2|) 



1600+ID 



99 /qq 



77° 



'ID as in IPROC = 300 + ID) 



1700+IQ 



QCD heavy quark production (IQ as above) 

After generation, IHPRO is subprocess (see sect. [4.6.2|) 



1800 QCD direct photon + jet production 

Table ^: Process codes. (Continues) 



36 



IPROC 


Process 




After generation, IHPRO is subprocess (see sect. 


4.6.51) 


1900+ID 


qq q'q'W+W- /Z^Z^ (ID as in IPROC = 300 + ID) 


2000 
2001-4 
2005-8 


t production via exchange (sum of 2001-2008) 
ub ^ di , db ^ ui , db ^ ui , ub ^ dt 
cb ^ si , sb ^ ci , sb ^ ct , cb st 


2100 
2110 
2120 
2150 
2160 
2170 


+ jet production 

+ jet production (Compton only: gq Wq) 
+ jet production (annihilation only: qq Wg) 
Z° + jet production 

Z° + jet production (Compton only: gq Zq) 
Z^ + jet production (annihilation only: qq Zg) 


2200 


QCD direct photon pair production 

After generation, IHPRO is subprocess (see sect. 


4.6.5|) 


2300+ID 


QCD SM Higgs + jet production (ID as in IPRO 
After generation, IHPRO is subprocess (see sect. 


C=30( 
4.6.10 


)+ID) 

) 


2400 
2450 


Mueller- Tang colour singlet exchange 
Quark scattering via photon exchange 


2500+ID 


gg/qq ttH^^ (ID as in IPR0C=300+ID) 


2600+ID 


qq' W^H^^ (ID as in IPR0C=300+ID) 


2700+ID 


qq Z^H^M (ID as in IPRGC=300+ID) 


2800 
2810 
2815 
2820 
2825 
2850 
2860 
2870 
2880 


W'^W" production in hadron-hadron collisions 

^0^0 production in hadron-hadron collisions (including photon terms) 

^0^0 production in hadron-hadron collisions {Z^ only) 

W^Z^ production in hadron-hadron collisions (including photon terms) 

W^Z^ production in hadron-hadron collisions (Z° only) 

hadron-hadron W~^W~X using MCONLO 

hadron-hadron Z^Z^X using MCONLO 

hadron-hadron W+Z°X using MCONLO 

hadron-hadron W'Z^X using MCONLO 


2900+IQ 
2910+IQ 


gg + qq QQZ^ for massless Q and Q (IQ=1. . .Q ioi Q = d . . .t) 
gg + qq ^ QQZ^ , for massive Q and Q (IQ=1. . . 6 for Q = d . . .t) 


3000-3999 
3000 
3010 
3020 
3030 


Minimal Supersymmetric Standard Model (MSSM) processes 
2-parton 2-sparticle processes (sum of those below) 
2-parton 2-sparton processes 
2-parton 2-gaugino processes 
2-parton — > 2-slepton processes 


3100+ISQ 


gg/qq qq*H^ (ISQ=IPR0C-3100 as from table 15) 


3200+ISQ 


gg/qq qq*h, H, A (lSQ=IPR0C-3200 as from table 


16) 



Table ^: Process codes. (Continues) 
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IPROC 

3310,3315 
3320,3325 
3335 

3350 

3355 
3360,3365 
3370,3375 



Process 
qq' 
qq' ■ 
qq' ■ 
qq- 
qq- 
qq- 
qq- 



^ W^h^, H^h^{al\ q, q' flavours - gauge bosons mediated only) 
y W^H\H^H^ (") 
(") 

W^H^ (Higgstrahlung and Higgs mediated) 
H^H^ (all q flavours — gauge boson mediated only) 
Z^h'>,A^h'> (") 



3410 
3420 
3430 
3450 



bg ^ b + ch. conj. 
bg ^ b + ch. conj. 
bg ^ b + ch. conj. 
bg ^ t H~ + ch. conj. 



3500 



bq bq'H + ch. conj. 



3610 
3620 
3630 



ml 99 



(light scalar Higgs) 
(heavy scalar Higgs) 
A^ (pseudoscalar Higgs) 



3710 
3720 



3900-99 



qq q'q'W+W- /Z^Z^ q'q'h^ 
qq q'q'W+W- /Z^Z^ q'q'H^ 



3810+IQ 99 + ^ QQhP (all q flavours in s-channel, IQ as usual for Q flavour) 

3820+IQ gg + qq-^ QQH'' C) 

3830+IQ gg + qq^QQA"^ {") 

3839 99 + ~* biH^ + ch. conjg. (all q flavours in s-channel) 

3840+IQ gg QQh^ (IQ as above) 

3850+IQ gg~^QQH^{") 

3860+IQ gg^QQA^i^") 

3869 gg biH^ + ch. conjg. 

3870+IQ qq —>■ QQh^ (all q flavours in s-channel, IQ as above) 

3880+IQ qq^QQH^i") 

3890+IQ qq-^QQA°{") 

3899 qq — * biH^ + ch. conjg. (all q flavours in s-channel) 



Reserved for other hadron-hadron MSSM processes 



4000-99 
4000 
4010 

4010+IN 
4020 

4020+IC 
4040 
4050 



R-parity violating supersymmetric processes via LQD 
single sparticle production, sum of 4010-4050 



Ujdk 
Ujdk 
Ujdk 
Ujdk 
Ujdk 



djdk 



Xwh ' djdk 
X'l^i, djdk - 
Xic^i, djdk ■ 
f+Z°, Ujdk 



X^Vi (all neutralinos) 
Xm^i (IN=neutralino mass state) 

X~^t (^11 charginos) 
> Xic^t (IC=chargino mass state) 
^ UiW+ and djdk ^ ifW' 



Ujdk ifh^/H^/A^, Ujdk UiH+ and djdk 
Table ^: Process codes. (Continues) 
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IPRDC I Process 

4060 Sum of 4070 and 4080 

4070 Ujdk — > uidm and djdk 

4080 Ujdk — ^ t^j^k ^i^k - 



in 



via LQD only 
J „^ , via LQD and LLE 



4100-99 

4100 

4110 
4110 +IN 

4120 
4120 +IC 

4130 

4140 

4150 

4160 



R-parity violating supersymmetric processes via UDD 
single sparticle production, sum of 4110-4150 
Uidj X^dk, djdk 



Uidj ■ 
Uidj - 
Uidj - 

1^1 dj J 

Uidj 



(all neutralinos) 
Xjj^tij(IN as above) 
^ X~di (all charginos) 
Xic^i (IC as above) 
gui 

> Uidj i*W 
dl^h^/H^/A^,djdk u*^hyHyA^,Uidj - 

djdk 



XiN^ki djdk 
X^Uk, djdk 
xtc'^k, djdk 
gdk, djdk — 
hlZ\ djdk - 



and djdk 



b*W~ 



u 



Uidr, 



djdk 



didm via UDD. 



4200-99 
4200 
4210 

4210+IQ 
4220 
4250 

4250+IL 
4260 
4270 
4271 
4272 
4273 



Graviton resonance production 
Sum of 4210, 4250 and 4270 
gg/qq ^ G ^ gg/qq (all partons) 
qq (IQ as above) 
99 

U (all leptons) 
£i (lL = l-6 for ^ 



99 /qq 
99 /qq 
99/qq 
99 /qq 
99/qq 
99/qq 
99/qq 
99/qq 
99/qq 



G 
G 
G 
G 
G 
G 
G 
G 
G 



e, z/^, etc.) 



77 

w^w-/z'zyHi^Hi^ 
w+w- 

-"SM-"SM 



5000 
5100+IQ 
5200+IQ 

5300 
5500 
5510,20 



Pointlike photon-hadron jet production (all flavours) 
Pointlike photon heavy flavour pair production (IQ as above) 
Pointlike photon heavy flavour single excitation (IQ as above) 
After generation, IHPRO is subprocess (see sect. [4.6.5|) 
Quark-photon Compton scattering 

Pointlike photon production of light {u, d, s) L=0 mesons 

S=0 mesons only, S=l mesons only 

After generation, IHPRO is subprocess (see sect. [4.6.5| ) 



6000 77 qq (all flavours) 

6000+ 1 Q 77 — s> gg (IQ as above 

6006+IL 77 ^ il (IL 

6010 77 ^ W+W 



1, 2, 3 for £ = e, /z, r) 



7000 



Baryon-number violating and other multi-l^^ processes 
Table p|: Process codes. (Continues) 
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IPROC 


Process 


7999 


generated by HERBVI package 


8000 


Minimum bias soft hadron-hadron event 


9000 
9000+IQ 

9010 
9010+IQ 


Deep inelastic lepton scattering (all neutral current) 
Deep inelastic lepton scattering (NC on flavour IQ) 
Deep inelastic lepton scattering (all charged current) 
Deep inelastic lepton scattering (CC on flavour IQ) 


9100 
9100+IQ 

9107 

9110 
9110+IP 

9130 
9140+IP 


Boson-gluon fusion in neutral current DIS (all flavours) 
Boson-gluon fusion in neutral current DIS (IQ as above) 
J/'>p + gluon production by boson-gluon fusion 
QCD Compton process in neutral current DIS (all flavours) 
QCD Compton process in NC DIS (IP=1-12 ioid-t,d-t) 
All C(as) NC processes (i.e. 9100+9110) 

Heavy quark production by charged- current boson-gluon fusion 
IP: 1 = sc, 2 = 6c, 3 = st, A = ht {-V ch. conj.) 


9500+ID 


W+W-/Z^Z^ H^^ in DIS (ID as in IPROC = 300 + ID) 


10000+IP 


as IPROC = IP but with soft underlying event 

(soft remnant fragmentation in lepton-hadron) suppressed 



Table ^: Process codes. 



4.2.1 Treatment of quark masses 

The extent to which quark mass effects are included in the hard process cross sec- 
tion is different in different processes. In many processes, they are always treated as 
massless: IPROC = 1300, 1800, 1900, 2100, 2300, 2400, 5300, 9000. In two processes 
they are all treated as massless except the top quark, for which the mass is correctly 
incorporated: 1400, 2000. In the case of massless pair production, only quark flavours 
that are kinematically allowed are produced. In all cases the event kinematics incor- 
porate the quark mass, even when it is not used to calculate the cross section. In 
two processes, quarks are always treated as massive: 500, 9100. Finally, in several 
processes, the behaviour is different depending on whether a speciflc quark flavour 
is requested, in which case its mass is included, or not, in which case all quarks are 
treated as massless. These are: IPROC = 100, 110, 120, QCD 2^2 scattering (1500 
vs. 1700+IQ), jets in direct photoproduction (5000 vs. 5100+IQ and 5200+IQ). In 
the case of IPROC = 2900,2910 one has the option of using the massless or massive 
matrix element. 

These differences can cause inconsistencies between different ways of generating 
the same process. The most noticeable example is in direct photoproduction, where 
one can use process 9130, which uses the exact 2^3 matrix element e + g ^ 
e + q + q, OT process 5000, which uses the Equivalent Photon Approximation (EPA) 
for e — s> e + 7 and the 2^2 matrix element for 7 + — > g + g. For typical HERA 
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kinematics, the EPA is valid to a few per cent, but the difference between the two 
processes is much larger, about 20% for PTMIN = 2 GeV. This is entirely due to the 
difference in quark mass treatments, as can be checked by comparing process 9130 
with processes 5100+IQ and 5200+IQ summed over IQ. 

4.2.2 Couplings 

The two-loop QCD couphng at scale Q is given by subroutine HWUALF with arguments 
lOPT = 1 and SCALE = Q. Threshold matching is performed at the quark mass scales 
Q — RMASS(i). Setting lOPT = initialises the coupling using the 5-flavour value 
= QCDLAM. Other values of I OPT are for internal use only. 

The electromagnetic coupling is given by HWUAEM(Q2) = eV(47r); it runs accord- 
ing to the prescription in ref. [57] with the hadronic term as given in ref. [58]. The 
parameter ALPHEM = HWUAEM(O), default value 0.0072993, provides the normahzation 
at the Thomson limit; it is used for all processes involving real photons. Photon emis- 
sion in parton showers and in the 'dead-zone' in e'^e~ processes can be enhanced by 
a factor of ALPFAC (default = 1). The normalised electric charges of the fundamental 
fermions are stored in the array QFCH(I), where I = 1-6 for the quarks d, u, s, c, b, t 
(e.g. QFCH(4) = 2./3.) and 11-16 for the leptons e, z/g, /i, z/^, r, z^^. 

The weak neutral current is taken to be of the form e(f j+a/75)7^, where the elec- 
tric charge is evaluated at a scale appropriate to the process. The arrays VFCH(I, J) 
and AFCHdjJ) store the couplings: I as before, J = 1 for the minimal Standard 
Model and 2 for possible Z' couphngs (only used if ZPRIME= . TRUE . ) . Note that 
universality is not assumed — couplings can be arbitrarily set separately for each 
fermion species. The default couphngs are given in terms of of SWEIN= sin^ 9w, 
default value 0.2319, as: 

Vf = (T3/2 — Q sin^ 9w) / (cos 9w si^n. 9w) , ctf = T3/(2cos6'vFsin6'iy) . 

The weak charged current is given in terms oi g — e/ sin 9w and the Cabbibo- 
Kobayashi-Maskawa mixing matrix, the elements squared of which are stored in 
VCKM(K,L), K = 1,2,3 for u,c,t, L = 1,2,3 for d,s,b. The variable SCABI = 
sin^ ^Cabibbo IS however also retained for the present. Note the Fermi constant Gpermi 
is eliminated from all cross sections. 

The overall scale for all cross sections, given in nanobarns, is set by GEV2NB = 
{hc/ef, default value 389379. 

We now give more detailed descriptions of the various subprocesses, concentrat- 
ing again on the new features since ref. [1] . 

4.3 Lepton-antilepton Standetrd Model processes 

Lepton beam polarization effects are included in e'^e" — > 2/3 jet production and the 
Bjorken process {ZH production). Incoming lepton and antilepton beam polariza- 
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tions are specified by setting the two 3-vectors EPOLN and PPOLN: component 3 is 
longitudinal and 1,2 transverse. 

Photon initial-state radiation (ISR) in e+e~ annihilation events is allowed. The 
parameter TMNISR sets the minimum s/s value (defauh = 10"^), ZMXISR sets the 
(arbitrary) separation between unresolved and resolved emission (default = 1 — 10^^). 
Setting ZMXISR = switches off photon ISR. 

Where processes are listed for they are available for e~^e~ and fi'^fi~ anni- 
hilation. Many of the processes listed for e+e" will also work for fi'^fj^~, but we have 
not been systematic in ensuring this. 



4.3.1 IPROC=100-127: hadron production 

A correction to hard gluon emission in e^e~ events has been added and is now the 
default process for IPR0C=100+IQ. The 0{as) matrix element is used to add events 
in the 'dead zone' of phase-space corresponding to a quark- ant iquark pair recoiling 
from a hard gluon [16]. Although this is asymptotically negligible, and cannot be 
produced within the shower itself, it has a significant effect at LEPl energies. As a 
result, the default parameters have been retuned, and show a marked improvement 
in agreement with e+e~ data for event shapes sensitive to three-jet configurations. 

The routine HWBDED implements this hard correction while HWBRAN has been 
modified to include the soft matrix-element corrections described in section B.2.3. 



When IPR0C=100+IQ, hard gluons emitted into the dead zone are assigned to 
the quark or antiquark shower and do not appear explicitly in the event record. 

The qqg process alone, generated according to the 0{as) qqg matrix element 
with a maximum thrust cutoff THMAX (defauh 0.9), is given by IPR0C=110+IQ. 

The uncorrected qq process has been retained for comparative purposes and is 
available as IPR0C=120+IQ. 

The fictional e+e^ processes e+e^ 9 + g{+g), IPR0C=107 and 127, is treated 
just like e"^e~ — > qq, summed over light quark flavours, for direct comparisons be- 
tween quark and gluon jets. 



4.3.2 IPR0C=15O— 250: lepton and electroweak boson production 

In IPR0C=150, only the s-channel process, mediated by a virtual photon or Z^, is 
included, so the final-state leptons must be different from the initial ones. 

The processes of W~^W~ and Z^Z^ pair production, IPR0C=200 and 250, are 
based on a program kindly supplied by Zoltan Kunszt, which fully includes decay 
correlations. The QCD 0{as) matrix element correction for hard gluon emission in 
hadronic W and Z decays has also been implemented in these processes, according 
to the method described in section p.2.3| . In contrast to IPR0C=100+IQ, any hard 
gluons emitted into the dead zone are shown explicitly in the event record. 
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4.3.3 IPROC=300-499: Higgs boson production 

HERWIG generates SM Higgs bosons in lepton-antilepton collisions through the 
Bjorken process Z^*) Z^*^Hly^ with one or both Z°'s off-shell (IPROC = 300 + ID) 
and W^W~ /Z^Z^ fusion (IPROC = 400 + ID). See section ^ for explanation of how 
the Higgs decay is controlled by the value of ID. 

4.3.4 IPROC=500— 559: two-photon/photon-boson processes 

In the e"'"e^ two-photon processes, IPROC = 500 + ID, ID = — 10 is the same as 
in Higgs processes for gg, ^^ and W^W^ . The Equivalent Photon Approximation 
(EPA) is used for the e ^ 67 vertices. The phase space is controlled by EMMIN and 
EMMAX for the two-photon centre-of-mass frame (CMF) mass, PTMIN and PTMAX for 
the transverse momentum of the CMF in the lab, and CTMAX for the cm. angle of 
the outgoing particles. The additional phase-space variable WHMIN sets the minimum 
allowed hadronic mass and affects photoproduction reactions (7-hadron and 7-7) and 
DIS. 

In photon-iy^ fusion, IPROC = 550 + ID, ID = — 9 is also the same as in Higgs 
processes, except that ID=1 or 2 both give the sum of du and ud etc. The EPA 
is used for the e — 67 vertex. The phase space is controlled by EMMIN and EMMAX 
only. The full 2 — > 3 matrix elements for 76 — > //V are used, so the cross section 
for real production is correctly included. In the case of 77 — > WW the decay 
correlations are not yet correctly included: the W^s currently decay isotropically. 

4.3.5 IPROC=600-656: four jet production 

Electron-positron annihilation to four jets is provided by IPR0C=600+IQ, where a 
non-zero value for IQ guarantees production of quark flavour IQ whilst IQ=0 cor- 
responds to the natural flavour mix. IPRGC=650+IQ is as above but without those 
terms in the matrix element which orient the event w.r.t. the lepton beam direction. 
The matrix elements are based on those of Ellis, Ross and Terrano [59] with orien- 
tation terms from Catani and Seymour [60]. The soft and collinear divergences are 
avoided by imposing a minimum ycnti Y4JT (default 0.01), on the initial four partons. 
The interjet distance ?/cut is calculated using either the Durham or JADE metrics, 
as selected by the logical variable DURHAM (default .TRUE.). In order to improve effi- 
ciency parameterizations of the volume of four-body phase space are used: these are 
accurate up to a few percent for Ucnt values less than 0.14. Note also that the phase 
space is for massless partons, as are the matrix elements, though a mass threshold 
cut is apphed. 

The argument of the strong coupling is set equal to EMSCA, the scale for the 
parton showers. This is taken to be the smallest of the ?/cut values times the cm. 
energy squared if FIX4JT=. FALSE, (default); otherwise the argument is fixed at Y4JT 
times the cm. energy squared. 
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The matrix elements for the qqgg and 
qqqq (same flavour quark) final states re- 
ceive contributions from two colour flows 
each. The actual process and colour flow 
generated is indicated by IHPRO as shown 
in table [TO. The meaning of 'c/f conn.' is 



IHPRO 


1+2+3+4 


c/f conn. 


91 


Q -\- Q -\- Q -\- Q 


3 14 2 


92 


q + q + 9 + 9 


4 12 3 


93 


q+q+q+q 


4 12 3 


94 


q+q+q+q 


2 14 3 


95 


q + q + q' + q' 


4 12 3 



Table 10: Four jet subprocesses. 



discussed in section [4.6.2| below. The treat- 
ment of the interference terms between the 
two colour flows is controlled by the array 

I0P4JT(1) for qqgg and I0P4JT(2) for qqqq (identical quark flavour): 



I0P4JTfl 



neglect 

1 extreme 3142 

2 extreme 4123 



I0P4JT(2) 



neglect 

1 extreme 4123 

2 extreme 2143 



In both instances the default value is 0. 

See ref. [61] for some applications and discussions of the new four-jet implemen- 
tation in e~^e~ annihilations. 

4.4 Lepton-antilepton supersymmetric processes (MSSM) 

The R-parity conserving lepton-antilepton SUSY processes have IPROC = 700 — 799 
and processes have IPROC = 800 — 899. Lepton beam polarization effects are 
included for all of the SUSY production processes. The processes have all been 
implemented in such a way as to allow either e~^e~ or as the initial state. As 

with the SM lepton-antilepton processes, ISR is allowed for the SUSY production 
processes. 

As, by probing the individual thresholds, it may be possible to study the produc- 
tion of a given sparticle pair in lepton-antilepton collisions, we have provided more 
control over the sparticles produced than for the hadron-hadron SUSY production 
processes described in section 

We remind the reader here that all SUSY particle data have to be read from an 



input file before event generation (see section 3.5.1 ). 



4.4.1 IPROC=700— 799: gaugino, slepton and/or squark production 

With IPROC = 700 one obtains the four processes IPROC = 710, 730, 740, 760 in 
the correct proportions. The matrix elements have been derived independently and 
the cross sections are in good agreement with those from SUSYGEN [62]. For all 
these processes the hard process scale EMSCA has been set to the centre-of-mass 
energy. 

The gaugino and sfermion mixing conventions of Haber and Kane [40] are adop- 
ted in all cases. In addition, the neutralino mixing matrix ZMIXSS is defined internally 
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in terms of the photino, zino and current eigenstate neutral higgsino components, 
instead of the bino, W^a-ino and higgsino components adopted for ZMXNSS, equivalent 
to Nij in [40]. 

IPROC = 710 gives lepton-antilepton— i>neutralino pair production. A number 
of additional IPROC codes have been provided to enable the user to produce given 
neutralino mass eigenstates. It should be noted that in order to provide a compact 
code for these processes there is more than one possible IPROC number for some 
processes. For example the final state X1X3 can be produced using either the IPROC 
codes 713 or 721. 

IPROC = 730 gives lepton-antilepton— >chargino pair production. As with the 
neutralino production there are codes to allow a given pair of charginos to be pro- 
duced. 

IPROC = 740 gives lepton-antilepton— i>slepton pair production. In these pro- 
cesses for the first two generations the left /right eigenstates are produced, while 
for the third generation the mass eigenstates are produced. In the processes pro- 
ducing given slepton pairs for f production processes the left eigenstate is replaced 
by the lighter mass eigenstate and the right eigenstate by the heavier mass eigen- 
state. 

IPROC = 760 gives lepton-antilepton— >squark pair production. As with the other 
processes additional codes are provided to allow the production of given qq* pairs. 
For stop and sbottom production the mass eigenstates are produced, while for the 
first two generations the left/right eigenstates are generated. As with the slepton 
production for the third generation when a given squark pair is requested the left 
eigenstate is replaced by the lighter mass eigenstate and the right eigenstate by the 
heavier mass eigenstate. 



4.5 Other lepton-antilepton non-Standard-Model processes 
4.5.1 IPROC=800-899: R-parity violating SUSY processes 

A range of possible production processes in lepton-antilepton collisions is included. 



Unlike the case of production in hadron-hadron collisions, section |4.8.1| , we have 
included processes for which there is either no s-channel resonance or the resonance 
is not kinematically accessible. 

All the possible single sparticle production mechanisms which occur via the first 
term in the superpotential given in [42] are included. This includes the process 
t^l~ —>■ 7z/ for which there is no s-channel resonance. As the cross section for this 
process diverges in the limit that the photon is collinear with the incoming lepton 
beams a cut on the pt of the outgoing particles pt > PTMIN has been imposed for 
this process (IPR0C=850). The ISR is switched off for this process as including it 
would lead to a double counting of the photon radiation. For this reason this process 
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is not included in the code IPR0C=800 which generates all the other single sparticle 
production mechanisms. 

We have also included some processes for the production of SM particles via 
s-channel sneutrino exchange. In addition to the s-channel sneutrino diagrams the t- 
channel sparticle exchanges, SM diagrams and all the interference terms are included. 
This uses a generahzation of the formulae of [63] . A cut pr > PTMIN is used in the 
process £~^£~ — > (IPR0C=870) to avoid the divergence as i — > in the Bhabha 
scattering cross section. 

Except where stated explicitly above, no PTMIN cut is applied. 

4.5.2 IPROC=900-999: MSSM Higgs production processes 

For MSSM Higgs pair production (lPR0C=955-965) a new subroutine has been intro- 
duced, HWHIHH, whereas for the other processes we make use of the implementation 
of their SM counterparts, which are based on the subroutines HWHIGW and HWHIGZ. 

4.5.3 IPR0C=1000-1199: SM and MSSM Higgs associated production 

SM and MSSM Higgs production in association with fermion pairs in lepton-lcpton 
collisions is available in version 6.5. These processes were introduced in [64] and their 
phenomenological relevance was discussed in [65]. 

Fermion masses are retained in the final state according to the HERWIG de- 
faults. The same values appear in the Yukawa couplings. Notice that in the case 
of charged Higgs boson production the Cabibbo-Kobayashi-Maskawa mixing matrix 
has been assumed to be diagonal. Furthermore, due to the rather different phase 
space distribution of the final state products, all processes can only be produced 
separately, not collectively. Initial- and final-state radiation (both QED and QCD) 
and bcamsstrahlung are included via the usual HERWIG algorithms. Finally notice 
that the use of the IPROC series 1000 and 1100 for i'^fr processes required some 
internal modification to HERWIG, which was originally designed to generate lep- 
tonic processes only for IPROC < 1000. Now we assume an i'^i" process whenever 
IPROC < 1300. These modifications have no imphcations for the traditional user, but 
may affect more knowledgeable ones who have edited previous versions of the main 
HERWIG code. 

4.6 Hadron-hadron Standard Model processes 
4.6.1 IPR0C=1300-1499: Drell-Yan processes 

The Drell-Yan code is extended to the production of all fermion pairs; 1300 gives all 
quark fiavours; 1300-l-IQ a specific quark fiavour, 1350 all leptons (including neutri- 
nos) 1350-l-IL a specific lepton fiavour. The s-channel component of the interference 
with like-fiavour qq scattering is included here. 
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IHPRO 


1 + 2 - 


3 + 4 


c/f conn. 


1 


q + q - 


q + q 


3 4 2 1 


2 


q + q - 


q + q 


4 3 12 


3 


q + q' - 


q + q' 


3 4 2 1 


4 


q + q - 


q' + q' 


2 4 13 


5 


q + q - 


q + q 


3 14 2 


6 


q + q - 


q + q 


2 4 13 


7 


q + q - 


g + g 


2 4 13 


8 


q + q - 


g + g 


2 3 4 1 


9 


q + q' - 


q + q' 


3 14 2 


10 


Q + g - 


q + g 


3 14 2 


11 


Q + g - 


Q + g 


3 4 2 1 


12 


q + q - 


q' + q' 


3 14 2 


13 


q + q - 


q + q 


2 4 13 


14 


q + q - 


q + q 


3 14 2 


15 


q + q - 


g + g 


3 14 2 


16 


q + q - 


g + g 


4 12 3 


17 


q + q' - 


^ q + q' 


2 4 13 


18 


q + q - 


q + q 


4 3 12 


19 


q + q - 


q + q 


3 4 2 1 


20 


q + q' - 


q + q' 


4 3 12 


21 


Q + g - 


q + g 


2 4 13 


22 


Q + g - 


q + g 


4 3 12 


23 


g + Q - 


g + q 


2 4 13 


24 


g + Q - 


g + q 


3 4 2 1 


25 


g + Q - 


g + q 


3 14 2 


26 


g + Q - 


g + q 


4 3 12 


27 


g + g - 


q + q 


2 4 13 


28 


g + g - 


q + q 


4 12 3 


29 


g + g - 


g + g 


4 12 3 


30 


g + g - 


g + g 


4 3 12 


31 


g + g - 


g + g 


2 4 13 



Table 11: QCD subpro cesses. 
The initial-state parton showers in Drell-Yan processes are matched to the ex- 



act 0{as) matrix-element result as discussed in section |3.2.3| . The routine HWBDYP 
implements the hard corrections whilst HWSBRN includes the soft corrections to the 
initial-state radiation. For further details see ref. [20]. 

4.6.2 IPROC=1500: QCD 2^2 processes 

At present only 2 — 2 subprocesses are implemented. They are classified in ta- 
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ble [TT|. Here and in other subprocess tables, 'c/f conn.' refers to the colour /flavour 

connections between the partons: H j k V means that the colour of parton 1 comes 

from parton i, that of 2 from j, etc. For antiquarks, which have no colour (only 

anticolour), the label shows instead to which parton the flavour is connected. For 

this colour /flavour labelling all partons are defined as outgoing. Thus, for example, 

process 10 has colour connections 3 14 2, corresponding to the colour flow diagram: 

1 3 



~2"^ ' 4 

When different colour flows are possible, they are listed as separate subprocesses. 
This separation is not exact but is normally a good approximation [6,7]. The sepa- 
ration is now performed using the improved method proposed in [31], as outlined in 



section 3.1. The sum of the colour flows is the exact lowest-order cross section. 



4.6.3 IPROC=1600— 1699: Higgs boson production by parton fusion 

IPROC = 1600 + ID gives the sum of gg and qq fusion. The lowest-order formulae 
that we have used can be found in [41]. The hard processes are implemented in the 
subroutine HWHIGS. 

4.6.4 IPROC=1700-1706: heavy quark production 

The separation of colour flows is now performed using the improved method proposed 



in [31], as outlined in section |3.1| . The classification of subprocesses according to 
IHPRO is as for IPR0C=1500. 

4.6.5 IPROC=1800: QCD direct photon plus jet production 

The relevant IHPRO codes are 41-47 in table 0. For future reference we also collect 
here the codes for other processes that involve outgoing direct photons or incoming 
pointlike photons (IPR0C=2200, 5000-5520). Note that the photon is colour /flavour- 
connected to itself. In the cases IHPR0=71-76, M represents an L = meson (see 
IPR0C=5500). 

4.6.6 IPROC=1900— 1999: Higgs boson production by weak boson fusion 

The qq q'''^(fWV q^'^^q'^'^^H^^ subprocesses, for VV = W^W^ , Z^Z^, summed 
over initial- and final-state quarks can be invoked by setting IPR0C=1900-t-ID, with 
ID used to identify the Higgs decay. 

The formulae used are well known in the literature and can be found e.g. in [22]. 
This process is administered by the subroutine HWHIGW, which also handles the similar 
cases initiated by e+e~ and e^p collisions. 
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4.6.7 IPROC=2000-2008: single top 
production 

The process of single top quark pro- 
duction by ly-boson exchange includes 
so far only those processes initiated 
by a 6-quark (or antiquark) and a first- 
or second-generation quark or anti- 
quark. Note that this requires 6-quarks 
to be available in the parton distribu- 
tion functions. 

4.6.8 IPROC=2100-2170: electro- 
weak boson plus jet production 

The electroweak boson decay corre- 
lations and width are now correctly 
included in these processes. 

4.6.9 IPROC=2200: direct photon 
pair production 

See section ^. 6. 5| for IHPRO codes (GI- 
GS). 



IHPRO 



41 
42 
43 
44 
45 
46 
47 



51 
52 
53 



61 
62 
63 



71 
72 
73 
74 
75 
76 



1 



4 



q + q 

1 + 9 

q + q 

Q + 9 

9 + (l 

9 + Q 

9 + 9 



9 + 1 
Q + 1 
9 + 1 
Q + 1 
Q + 1 
Q + 1 
9 + 1 



1 + Q 
1 + Q 
1 + 9 



9 + Q 
9 + Q 
Q + Q 



Q + Q 
Q + Q 
9 + 9 



1 + 1 
1 + 1 
1 + 1 



7 + g ^ M{S = 0) + g' 
^ + q^ M{S = l)L + q' 
7 + g M{S = l)T + q' 
7 + g ^ M{S = 0) + g' 
7 + g ^ M{S = 1)l + q' 
-f + q^ M{S = 1)t + q' 



c/f conn. 



2 3 14 

3 12 4 
3 12 4 
2 3 14 

2 3 14 

3 12 4 
2 3 14 



14 2 3 

13 4 2 

14 2 3 



2 13 4 
2 13 4 
2 13 4 



14 3 2 

14 3 2 

14 3 2 

14 3 2 

14 3 2 

14 3 2 



Table 12: Direct photon subprocesses. 



4.6.10 IPROC=2300-2399: Higgs boson plus jet production 

High transverse momentum, scalar Higgs production via one-loop diagrams, in asso- 
ciation with a jet, is available as IPR0C=2300, within the SM. Only the top quark 
is included in the loops with lAPHIG controlling the approximation used: IAPHIG=0 
gives the zero top mass limit, 1 (default) the exact result, 2 the infinite top mass 
limit. The various subprocesses are illustrated in table ff^. 

Note that the Higgs boson is co- 
lour/flavour connected to itself. 

The relevant routines HWHGJl , 
HWHGJA/B/C/D, HWUCI2 and HWULI2 
use (non-standard Fortran 77) DOUBLE 
COMPLEX variables which may not be 
accepted by some compilers and are 
called COMPLEX* 16 by others. Users 
can change to COMPLEX variables, but 
this involves a risk of rounding errors 
spoiling numerical cancellations. 



IHPRO 


1+2 ^ 3+4 


c/f conn. 


81 

82 
83 
84 
85 
86 
87 


Q + Q 9 + H^M 
Q + 9 Q + 
Q + Q 9 + ^SM 
Q + 9 Q + ^SM 
9 + Q q + H^M 
9 + q q + H^M 
9 + 9 9 + H^M 


2 3 14 

3 12 4 
3 12 4 
2 3 14 

2 3 14 

3 12 4 
2 3 14 



Table 13: Higgs plus jet subprocesses. 
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4.6.11 IPROC=2400-2450: colour singlet exchange 

IPR0C=2400: Two-to-two parton scattering via exchange of a colour singlet, Mueller- 
Tang pomeron [66]. The fixed as and uq are given by ASFIXD (default 0.25) and 
OMEGAO (0.3) respectively. 

IPR0C=2450: Photon exchange, for like-flavour qq pairs including the t-channel com- 
ponent of the interference with qq —>■ qq via an s-channel photon or Z^. 

4.6.12 IPROC=2500— 2599: Higgs boson plus top quark pair production 

The SM 2 — > 3 Higgs production subprocesses of the type gg tiH^y^ and qq — > 
ttH^y^, for any flavour of initial state quarks q, are new to HERWIG version 6 and 
are handled by the subroutines HWHIGQ and HWH2QH. They are invoked by setting 
IPR0C=2500+ID (both gg and qq), with ID administering the Higgs decay channels 
as in IPR0C=300+ID. The initial state quark flavours q are always summed over. 
Notice that, given the size of the Yukawa couplings of the Higgs boson to quarks, 
in practice only the associated production with top quarks is of phenomenological 
relevance in the SM. The matrix elements used in the implementation can be found 



in ref. [22]. The treatment of the Higgs width here is as described in section |3!4 . 

4.6.13 IPROC=2600-2799: Higgs plus weak boson production 

The associated production of SM Higgs scalars with (lPR0C=2600-2699) and 
Z° (lPR0C=2700-2799) gauge bosons initiated by quark- ant iquark fusion via the 
2 — > 2 processes qq W^H^^ and qq Z^H^j^ is now available. A summation 
is as usual intended on the incoming quarks. The formulae given in [22] are used. 
The production processes are generated by the two new subroutines HWHIGV and 
HWH2VH whereas the Higgs decays are administered through the ID increment, as in 
IPR0C=300+ID. Again, the treatment of the Higgs boson width here is as discussed 



in section 



4.6.14 IPROC=2800-2825: Gauge boson pair production 

In version 6.3, the code already included in HERWIG for e^e^ —>■ WW/ZZ [67] was 
adapted for hadron-hadron collisions, including photons and the photon/Z interfer- 
ence for the resonant diagrams. 

All of these processes use a cut EMMIN (default value 20 GeV) on the mass of 
the gauge bosons produced. The cut PTMIN (default 10 GeV) on the transverse 
momentum of the bosons is also used. Both these cuts should not be taken to zero 
simultaneously if photon terms are included. The phase space for these processes 
contains a number of peaks and it was therefore necessary to use an adaptive multi- 
channel phase-space integration method which is described below. 

In versions 6.4 and above, the hard process scale EMSCA for this process has been 
changed to the parton-parton centre-of-mass energy from the average of the produced 
boson masses, which was used previously. 
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4.6.15 IPROC=2850-2880: Gauge boson pairs using MC@NLO 

These codes activate the interface to the program MCONLO for diboson production 
at next-to-leading order (see sect. [9l^ ). 



4.6.16 IPROC=2900-2916: QQZ production 

The matrix elements of [68] were used for the massless case and an independent 
calculation, using the approach of [69] , which was checked both for gauge invariance 
and against the massless case for the massive result. 

In both cases the decay of the Z is fully included and is selected using MODBOS. 
PTMIN controls the minimum transverse momentum of the outgoing quarks. As with 
gauge boson pair production it was necessary to use an optimized multi-channel 
phase-space integrator which is described below. 

The phase space for both gauge boson pair production and QQZ is complicated, 
as these processes are both treated as 2 — > 4 processes. In order to obtain a reasonable 
efficiency it was necessary to adopt a multi-channel approach based on that described 
in [70,71]. 

In each case a number of different channels are included which attempt to map 
the phase space for the different processes. The default weights for these different 
channels have been chosen to optimize the efficiency for the Tevatron and LHC; a 
choice of which to use is made based on the beam energy. However, the choice is 
affected by the phase space cuts applied. Therefore if these are significantly altered 
the weights for the different channels need re-optimizing. 

This is controlled by the new variable OPTM (defauh .FALSE.). If OPTM=.TRUE., 
before performing the initial search for the maximum weight HERWIG will attempt 
to optimize the efficiency using the procedure suggested in [71]. 

This is done by generating lOPSTP (default 10) iterations of lOPSH (default 1000) 
events. The choice of lOPSTP and lOPSH is a compromise between run time and accu- 
racy. The value of lOPSH should not be significantly reduced because the procedure 
attempts to minimize the error on the Monte Carlo evaluation of the cross section and 
if lOPSH is small the error on the error can be significant. If you need to re-optimize 
the weights we would recommend a long run just to optimize the weights which can 
then be used in all the runs to generate events. The new subroutine HWIPHS was 
added to initialise the phase space. 

4.7 Hadron-hadron supersymmetric processes (MSSM) 

The R-parity conserving SUSY processes occupy the IPROC entries of the series 3000 
(sparticle processes) and 3300-3600 (Higgs boson production), while processes 
have IPROC = 4000 - 4199. 

As with the lepton-antilepton SUSY processes the SUSY particle data must be 



read in from an input file before event generation as described in section p.5.1|. In 
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particular, unlike those of the SM Higgs boson H^^, the widths and decay modes of 
the SUSY Higgs bosons are not computed by HERWIG. 

A large number of final states involving the production of both neutral and 
charged Higgs bosons of the Minimal Supersymmetric Standard Model (MSSM) 
have been made available in version 6.3. They all proceed via 2^3 body hard 
scattering subprocesses. They are listed below, with corresponding process numbers 
(IQ and ISQ are as detailed in the following subsections). Further details of their 
implementation can be found in [22] 

The new subroutines introduced to administer the following processes are HWHIBQ 
and HWH2BH for IPR0C=3500, plus HWHISQ and HWH2SH for IPR0C=3100, 3200. In 
addition, HWHIGQ has been modified to accommodate the IPR0C=3800 series. 

4.7.1 IPROC=3000— 3030: sparton, gaugino and/or slepton production 

With IPR0C=3000 one obtains the three following processes, IPR0C=3010, 3020, 
3030, in the correct proportions. The variable IHPRO gives the subprocess actually 
generated (table 0). 

The matrix elements have been derived independently and the cross sections are 
in good agreement with those from ISAJET [39]. 

The hard process scale EMSCA has to be chosen globally for all sparton processes, 
e.g. for the argument of the QCD coupling. This is done using the kinematics ap- 
propriate to production of the lightest supersymmetric particle (LSP) and 



EMSCA 



2st'u' 



+ t'2 + u'2 

with i' = i — m^, u' = u — where m is the LSP mass. 

IPROC = 3010 gives 2-parton — > 2-sparton processes. All QCD sparton, i.e. 
squark and gluino, pair production processes are implemented. The matrix elements 
and the scheme for separating different colour flow parts are as given in [31]. 

IPROC = 3020 gives 2-parton 2-gaugino or gaugino-|-sparton processes. All 
gaugino, i.e. chargino and neutralino, pair production processes and gaugino-sparton 
associated production processes are implemented. 

The gaugino and sfermion mixing conventions of Haber and Kane [40] are used 
as described in section 

The various subprocesses, which include the 1 2 and charge conjugate reac- 



tions omitted below for brevity, are shown in table |I4. 

Note that the gauginos connect to themselves. The indices a, 6, z, j label gauginos 
in the order of increasing mass and take values a,b = 1 — 2 and i,j = 1 — 4. 

Gaugino mixing matrices are implemented in all subprocesses. The associated 
production subprocesses IHPR0=3026, 3027 include stop and sbottom left-right mix- 
ings. CKM mixing is implemented in subprocesses IHPR0=3023, 3025, 3027 but 
neglected in subprocess IHPR0=3021. 
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IPRGC = 3030 gives 2-parton 2- 
slepton processes. All Drell-Yan slep- 
ton production processes are imple- 
mented. The formulae agree with 
those of refs. [72, 73] in the limit of 
no stau left-right mixing. 

4.7.2 IPROC=3110-3178: charged 
Higgs plus squark pair produc- 
tion 



IHPRO 


1+2 ^ 3+4 


c/f conn. 


3021 
3022 
3023 
3024 
3025 
3026 
3027 


Q + Q xt + Xb 
Q + Q X° + X° 
a -\- a' — > y ^ + 

Q + Q X°i+9 

Q + q' xt + 9 

9 + Q X^ + Q 
9 + Q xt + Q' 


2 13 4 
2 13 4 
2 13 4 
2 4 3 1 
2 4 3 1 
2 4 3 1 
2 4 3 1 



Table 14: SUSY subprocesses. 



The production of charged Higgs bosons 

of the MSSM in association with squark pairs, of bottom and top flavours only, is im- 
plemented via the 3100 series of IPROC numbers (table |T5|). Their phenomenological 
relevance has been discussed in [74]. 

One should add 30(60) to IPROC for 
gg{qq)-on\j initiated processes. 



4.7.3 IPR0C=321O-3298: neutral 
Higgs plus squark pair production 

The production of neutral Higgs bosons 
of the MSSM in association with squark 
pairs, of bottom and top flavours only, 
is implemented via the 3200 series of 
IPROC numbers (table ^). Their phe- 
nomenological relevance has been dis- 
cussed in [74, 75]. 



IPROC 


partons - 


-> spartons 


Higgs 


3110 


99 + QQ - 


- Q^Q* 




3111 


99 + QQ - 




H+ 


3112 


99 + QQ - 


hi; 


H+ 


3113 


99 + QQ - 


hil 


H+ 


3114 


99 + QQ - 


U; 


H+ 


3115 


99 + QQ - 


^ ilk 


H 


3116 


99 + QQ - 


ilk 


H 


3117 


99 + QQ - 


i2k 


H 


3118 


99 + QQ - 


i2k 


H 



Table 15: Processes for IPR0C=3110-3178. 



One should add 30(60) to IPROC for gg{qq)-orAj initiated processes. 



4.7.4 IPROC=3310— 3375: Higgs-Higgs and Higgs-gauge boson pair 
production 

The production of Higgs-Higgs and Higgs-gauge boson pairs of the MSSM is im- 
plemented at tree level. We include gauge boson mediated contributions but not 
Higgstrahlung from the initial state, the only exception being W^H"^ production 
(IPR0C=3350), which does include diagrams where the Higgs boson couples to the 
initial partons as well as those mediated by neutral Higgs states [76]. Some of these 
processes are the MSSM equivalent of IPR0C=2600 and 2700 described earlier for the 
case of the SM. Here, the cases IPR0C=3310(3320) correspond to iy±/i°(iy±if°) and 
IPR0C=3360(3370) to Z^h^{Z^H^) final states. (No similar production can occur 
at leading order.) The array ENHANC is used to implement the MSSM couplings of 
Higgs scalars to gauge bosons. 
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IPROC 


partons - 


-> spartons 


Higgs 


3210(3220) [3230] 


gg + QQ - 




h{H)[A\ 


3211(3221) [3231] 


gg + qq - 




hiH)[A] 


3212(3222) [3232] 


gg + QQ - 


^ hb* 


h{H)[A] 


3213(3223) [3233] 


gg + QQ - 


^ kk 


h{H)[A] 


3214(3224) [3234] 


gg + QQ - 


hk 


h{H)[A] 


3215(3225) [3235] 


gg + QQ - 


ilil 


hiH)[A\ 


3216(3226) [3236] 


gg + QQ - 


tit*2 


h{H)[A] 


3217(3227) [3237] 


gg + qq - 


i2il 


h{H)[A] 


3218(3228) [3238] 


gg + qq - 




h{H)[A] 



Table 16: Processes for IPR0C=3210-3298. 



4.7.5 IPROC=3410— 3450: Higgs boson plus heavy quark production 

We have included so far only those processes initiated by a 6-quark (or antiquark) and 
a gluon. Note that this requires 6-quarks to be available in the parton distribution 
functions. 

4.7.6 IPROC=3500: charged Higgs boson from 6g-initiated processes 

This process is relevant for charged Higgs scalar production at large tan /? values, 
see [77]. 

4.7.7 IPROC=3610— 3630: neutral Higgs production by parton fusion 

These processes are the MSSM analogues of the SM processes 1600 etc. Recall 
however that the MSSM Higgs decay modes are controlled by the SUSY input data 



(see section p.5.1| ) and not by the value of IPROC. The subroutines HWHIGS and HWHIGT 
have been modified to take account of squark loop contributions and parity-violating 
Higgs-fermion couplings in the MSSM case. 

4.7.8 IPROC=3710— 3720: neutral Higgs production via weak boson fusion 

These processes are the MSSM counterparts of the SM process of weak vector-vector 
fusion in hadronic collisions (IPR0C=1900). They are computed using the same 
subroutines and setting the couplings to sin(/3 — a) for = (IPR0C=3710) 

and to cos(^ - a) for <l>° = (IPR0C=3720), respectively, where V = W^,Z^. 
(There is no A^VV coupling at tree level.) 

These reactions are two of the major direct production channels of neutral CP- 
even Higgs bosons of the MSSM at hadron colliders, such as the LHC (see e.g. [78]). 

4.7.9 IPROC=3811— 3899: Higgs boson plus heavy quark pair production 

For the case of neutral Higgs states, IPR0C=3810-|-IQ corresponds to production, 
IPR0C=3820+IQ to and IPR0C=3830+IQ to (For the last case, the vari- 
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able PARITY is set to —1.) Note also the production of charged Higgs states, via 
IPR0C=3839, 3869 and 3899, in association with pairs of top-bottom quarks. 

In the usage of the IPROC numbers corresponding to neutral Higgs states, when 
6-quarks are involved in gg-iusion modes (IPROC (+30) =3845, 3855 or 3865), the user 
should take care to avoid double-counting the chosen process with the corresponding 
2 ^ 1 and 2 ^ 2 cases of IPROC=3610-3630 and I PR0C=341 0-3430 initiated by 
quark- antiquark annihilation, i.e. bb Higgs, and (anti)quark-gluon scattering, i.e. 
bg ^ b Higgs, respectively: see [79]. Similar arguments hold for the charged Higgs 
states, as the (jf^f-induced process (IPR0C(-|-30)=3869) is an alternative implementa- 
tion of IPR0C=3450 [80]. 

The associate production of neutral Higgs bosons (both CP-even and CP-odd) 
of the MSSM with heavy QQ pairs [Q — h and t) is of extreme phenomenological 
relevance as a discovery mode of Higgs scalars, both at the Tevatron (Run 2) and 
the LHC (see e.g. [78,81]), as is the case of the charged Higgs channel [80,82]. 

4.8 Other hadron-hadron non-Standard-Model processes 

4.8.1 IPR0C=4000-4199: R-parity violating SUSY processes 

We include all the possible production processes of resonant sleptons and squarks 
in hadron-hadron collisions, for arbitrary numbers of non-zero couplings. These 
processes are implemented as two-to-two processes, i.e. with the decay of the resonant 
particle included. This allows us to include the t-channel diagrams where these 
occur. However we have not implemented those processes which can only occur via 
a t-channel diagram, or where the resonance will never be accessible. So for example 
while we include the process Uidj — > b\Z^ , which can occur via a resonant h^i we do 
not include the process Uidj — > h^Z^ , which cannot occur via a resonant diagram. In 
all cases both the processes listed and their charge conjugates are included. 

The scale choice is rather than the conventional transverse mass, due to the 
large number of different processes which must be calculated. The colour connection 
structure of these processes and their matrix elements can be found in ref. [42]. 

4.8.2 IPR0C=4200-4299 : graviton resonance production 

In some models with extra dimensions, Kaluza-Klein excitations of the graviton can 
be produced with significant cross sections at TeV-scale colliders. When the scale of 
the extra dimensions is not large, the excitations are manifest as discrete resonances. 

The production of a resonant excitation of the graviton is implemented as a 
2 — > 2 process including the decay of the resonance. The process is treated in 
a model-independent way, assuming only that there is a universal coupling of the 
graviton resonance to the SM fields. The effective lagrangian is given by 
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where /i^'^ is the spin-2 field and T^^, is the energy-momentum tensor of the SM fields. 
Although in these models there are usually many resonances, we have implemented 
only one. Others can be studied by making appropriate changes in the parameters. 
Graviton resonance production is described in more detail in [83]. 

The process is controlled by the coupling GRVLAM = A^, with dimensions of mass 
and default value 10 TeV, and by the mass EMGRV (default 1 TeV) and width GAMGRV 
of the resonance. If the width is set to zero (the default), the subroutine HWHGRV 
which calculates the cross section also calculates the width. 

The parton-level cross section for this process is non-unitary and is proportional 
to s/EMGRV^ at high energies. The fall-off of the parton distribution functions is not 
sufficient to suppress this bad high energy behaviour. Hence the parameters EMMIN 
and EMMAX controlling the minimum and maximum values of must be set. The 
default is to set these to 90% and 110% of the graviton resonance mass respectively. 
If the width of the resonance is more than a few percent of its mass then these limits 
should be reset. 

After event generation, IHPRO is set to 50 for qq initiated subprocesses and to 51 
for gg initiated subprocesses. 



4.9 Photon-hadron and photon-photon processes 

4.9.1 IPROC=5000— 5520: pointlike photon-hadron processes 

Pointlike photon-hadron scattering to produce QCD jets is available as IPROC = 
5000 — 5206. This is suitable for fixed-target photoproduction, provided events are 
generated in a frame in which the target has high momentum, and then boosted back 
to the lab. This is done if USECMF = .TRUE., the default, in which case the frame 
for event generation is the beam-target cm. frame. IPROC = 5100 -|- IQ gives fiavour 
IQ pair production, 7 + — > QQ, and IPROC = 5200 + IQ gives fiavour IQ single 
excitation, 7 + Q ^ 9 + Q, both including quark masses. IPROC = 5000 gives a 
sum over all processes and fiavours, 5100 and 5200, with massless quark kinematics. 
In all cases, after event generation the code IHPRO is set to 51, 52 or 53 according 



to the hard subprocess, as specified in section |4.6.5| . IPROC = 5300 gives Compton 
scattering, 7 + g — 7 + g. 

The direct, higher twist, production of light {u, d, s) L=0 mesons by point-like 
photons is also available: IPROC = 5500 for all spin =0 and 1 mesons; 5510 for only 
S=0 mesons; and 5520 for only S=l mesons. The vector mesons are produced with 
transverse or longitudinal polarization and decayed accordingly. The corresponding 
IHPRO codes (71-76) are also listed in section [4.6.5| . 

All these processes are available with lepton as well as hadron beams, using 
the Equivalent Photon Approximation. The phase-space variable WHMIN sets the 
minimum allowed hadronic mass and affects photoproduction reactions (7-hadron 
and 7-7) and DIS. In lepton-hadron DIS it is largely irrelevant since there is already 
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a cut on Bjorken y which at fixed s is almost the same, but for lepton-gamma DIS 
it makes a big difference. Direct 7 + 7* — g + g is included in the hard correction 
for lepton-gamma DIS. 

4.9.2 IPROC=6000— 6010: pointlike photon-photon processes 

Direct 77 — > charged particle pairs has been implemented with IPROC = 6000 + IQ: 
if IQ = 1 — 6 then only quark flavour IQ is produced, if IQ = 7, 8 or 9 then only 
lepton flavour e, /i or r is produced and if IQ = 10 then only W pairs are produced: 
in these cases particle mass effects are included. If IQ = 0, the natural mix of quark 
pairs is produced using massless matrix elements but including a mass threshold cut. 
The range of allowed transverse momenta is controlled by PTMIN and PTMAX as usual. 

4.10 Baryon number violating processes 

4.10.1 IPROC=7000-7999: generated by the HERB VI package 



See section 9.8 for details 



4.11 Minimum bias soft hadron-hadron collisions 

4.11.1 IPROC=8000: minimum bias soft hadron-hadron event 

Non-diffractive, soft hadronic, minimum bias events (IPR0C=8000) can be generated 
for the following combinations of beam and target: p,p, vr^, K^, e^, /i^, 7 on target 
p (or vice versa); p,p on target n (or vice versa); or 7 on 7. The event weight is 
the estimated cross section based on the parameterizations of Donnachie and Land- 
shoff [84]. The non-diffractive cross section is assumed to be 70% of the total. For 
lepton beams a photon is flrst generated using the Effective Photon Approximation 
(see section ^?T|) and then the on-shell photon cross section is used. See section p. 7. 3 
for discussion of the model used and the relevant parameters. 

4.12 Deep inelastic lepton scattering 

Deep inelastic (DIS) processes are broadly divided into those that start at 0(«s) 
(IPR0C=90**) and those like heavy quark and dijet production which start at 0{a^ 
(IPR0C=91**). Note that the DIS C(as) jet production processes, IPR0C=92**, have 
been withdrawn since they are subsumed within IPR0C=91**. 

The default limits on in DIS processes (Q2MIN, Q2MAX) have been set very 
small/large (0, 10^° GeV) and are reset to the kinematic limits unless changed by the 
user. This means the default Q2MIN is not suitable for simple neutral current DIS 
(IPR0C=9000 etc), but is appropriate for jet and heavy quark photoproduction. The 
range of the Bjorken-y variable {y = Q'^/xs) can be limited by YBMIN and YBMAX. 

The kinematic reconstruction of DIS processes can now take place in the Breit 
frame, if BREIT= . TRUE . (the default value). Previous versions used the lab frame. 
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Although the reconstruction is fully invariant under Lorentz boosts along the in- 
coming hadron's direction, it is not under transverse boosts, so there should be 
some difference between the two frames. The boost is not performed for very small 
Q^{< 10^^) to avoid numerical instabilities, but the two frames are in any case 
equivalent for such small values. 

The phase space for boson-gluon fusion is controlled by the parameters EMMIN, 



EMMAX (see section [4.3.4|) . The default values (0, ^/s) correspond to the behaviour of 
version 5.1. 

Lepton beam polarization effects are included in all DIS processes apart from J / ip 
production. The polarization is specified as in lepton-antilepton processes, i.e. by set- 
ting the 3- vectors EPOLN and PPOLN: component 3 is longitudinal and 1,2 transverse. 
Transverse only occurs in e'^e~ routines; recall that two transverse 'measurements' 
are needed to see an effect so it should not arise elsewhere. Note that in DIS processes 
one has to set either EPOLN if it is a lepton or (exclusive) PPOLN if an antilepton. 

All the DIS processes IPR0C=9000-9599 are available in e+e" as well as lepton- 
hadron collisions. The program generates a photon from the second beam (only) in 
the Equivalent Photon Approximation and the default is to use Drees-Grassie [53,54] 
structure functions for DIS on the photon. 

The parameter WHMIN sets the minimum allowed hadronic mass in DIS. In lepton- 
hadron DIS it is largely irrelevant since there is already a cut on Bjorken y which at 
fixed s is almost the same but for lepton-gamma DIS it makes a big difference. 

In addition to the processes listed here, a full simulation of QCD instanton-in- 
duced processes in DIS [85] is available through an interface to the program QCDINS 
[86] . For details, see the web page pittp : //www . desy . de/~tOOf ri/qcdins/qcdins . html| 



4.12.1 IPROC=9000-9006: neutral current 

Matrix-element corrections to DIS are available, following the general method de- 



scribed in section |3.2.3| . The hard correction is implemented in HWBDIS and the soft 
correction is included in routines HWBRAN and HWSBRN for the final- and initial-state 
radiation respectively. 

4.12.2 IPROC=9010-9016: charged current 

These are the charged current processes corresponding to those above, with the same 
treatment of hard and soft matrix element corrections. 

4.12.3 IPROC=9100-9130: 0{as) neutral current processes 

The photoproduction processes have been extended from the original heavy quark 
production program, to include all quark pair production (IPR0C=9100-9106) and 
QCD Compton (IPR0C=91 10-9122), as well as the sum of the two (IPR0C=9130). 
The possible flavours for processes 9100, 9110 and 9130 are limited by the input 
parameters IFLMIN and IFLMAX (defaults are 1 and 3, i.e. only u,d,s flavours). 
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A sign error has been corrected that led to the incorrect sign for the lepton-jet 
azimuthal correlation in QCD Compton processes in versions before 5.7. 

J /if) production (IPR0C=9107) now uses the Equivalent Photon Approximation 
instead of Weizsacker- Williams, with the same phase-space cuts as hadronic processes 
with lepton beams (see section |4?1|). 

The argument of the running coupling is controlled by the parameter BGSHAT — 
see section 

4.12.4 IPROC=9140— 9144: charged-current heavy quark production 

At present only W^g fusion processes are implemented. 



4.12.5 IPROC=9500— 9599: Higgs production by weak boson fusion 

The process of W^W" / Z^Z^ fusion into the SM Higgs boson is also available in e^p 
collisions, as IPROC = 9500+ID. For details of the implementation, see section [3^ and 
the corresponding processes initiated by lepton-lepton and hadron-hadron collisions, 
IPR0C=400+ID and IPR0C=1900+ID, respectively. 



4.13 Including new subprocesses 

The procedure for including further subprocesses remains substantially as described 
in ref. [1] but is repeated here for completeness. However, we now recommend that 
users make use of the Les Houches interface, see Section 

The parton and hard subprocess 4-momenta, masses and identity codes need 
to be entered in COMMON/HEPEVT/ with the appropriate status codes ISTHEP(I) = 
110-114 to tell the program which is which (see the table in section ^.3.1|) . The 



colour /flavour structure should be specified by the second mother and daughter 



pointers as explained in section [4.6.2 



The HERWIG identity codes IDHW(I) in COMMON/HWEVNT/ also need to be set cor- 
rectly. The IDHW codes can be listed in a run with IPRINT = 2: the most important 
are the quarks 1-6 (as IDHEP), antiquarks 7-12, gluon 13, overall centre-of-mass 14, 
hard centre-of-mass 15, soft centre-of-mass 16, photon 59, leptons 121-126, antilep- 
tons 127-132. 

The utility subroutine HWUIDT(IOPT, IPDG, IHWG,NAME) is provided to trans- 
late between Particle Data Group code IPDG, HERWIG code IHWG, and HERWIG 
CHARACTER*8 NAME, with lOPT = 1,2,3 depending on which of IPDG, IHWG and NAME 
is the input argument. 

Consider for example the process of virtual photon-gluon fusion to make 6 + 6 
in proton-electron collisions (in fact this process is included as IPROC = 9105). We 
assume the user provides a subroutine to generate the momenta PHEP for the hard 
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IHEP 


Entry 


ISTHEP 


IDHEP 


JMOHEP 


JDAHEP 


IDHW 


1 


e beam 


101 


11 








121 


2 


p beam 


102 


2212 








73 


3 


ep cm. 


103 











14 


4 


e in 


111 


11 


6 7 


7 


121 


5 


gluon 


112 


21 


6 9 


8 


13 


6 


hard cm 


110 





4 5 


7 9 


15 


7 


e out 


113 


11 


6 4 


4 


121 


8 


b 


114 


5 


6 5 


9 


5 


9 


b 


114 


-5 


6 8 


5 


11 



Table 17: Event record entries for eg ebb. 
subprocess e + g ebb. The colour structure is 



9 b 

Thus the momenta generated, together with those of the initial beams and the 



overall centre of mass, could be entered in the sequence shown in table |T7. 

Note that if there are more than two outgoing partons, the first has status 113 
and all the others 114. Each parton has JMOHEP(l, I) = 6 to indicate the location 
of the hard cm. for this subprocess, while JM0HEP(2, 1) gives the location of the 
colour mother (treating the incoming gluon as outgoing) or the connected electron. 
JDAHEP(1, 1) will be set by the jet generator HWBGEN, while JDAHEP(2, 1) points to the 
anticolour mother (or connected electron). Finally the HERWIG identifiers IDHW(I) 
could be set to the indicated values by means of the translation subroutine HWUIDT 
as follows: 

CHARACTER*8 NAME 



NHEP=9 
IDHEP(1)=11 
IDHEP (2) =2212 



IDHEP (9) =-5 
DO 10 I=1,NHEP 
1 CALL HWUIDT ( 1 , IDHEP ( I ) , IDHW ( I ) , NAME) 
IDHW(6)=15 

The last statement is needed because IDPDG(l) = returns IDHW(I) = 14. If subrou- 
tine HWBGEN is now called, it will find the coloured partons and generate QCD jets 
from them. Subsequent calls to HWCFOR etc. can then be used to form clusters and 
hadronize them. 
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If the hard subprocess routine is called from HWEPRO, like those already provided, 
it must have two options controlled by the logical variable GENEV in COMMON/HWHARD/. 
For GENEV = .FALSE., an event weight (normally the cross section in nanobarns) 
is generated and stored as EVWGT in COMMON/HWEVNT/. If this weight is accepted 
by HWEPRO, the subroutine is called a second time with GENEV = .TRUE, and the 
corresponding event data should then be generated and stored as explained above. 
On certain computers it will be necessary to SAVE those variables that determine 
event characteristics between the two subroutine calls. 

The parameter NMXJET sets the maximum number of outgoing partons in a hard 
subprocess (default 200). 



5. Parameters 



The quantities that may be regarded as adjustable parameters are indicated in ta- 



ble |T8|. Notes on parameters are given below. 



Name 


Description 


Default 


QCDLAM 


Aqcd (see below) 


0.18 


RMASS(l) 


Down quark mass 


0.32 


RMASS(2) 


Up quark mass 


0.32 


RMASS(3) 


Strange quark mass 


0.50 


RMASS(4) 


Charmed quark mass 


1.55 


RMASS(5) 


Bottom quark mass 


4.95 


RMASS(6) 


Top quark mass 


174.3 


RMASS(13) 


Gluon effective mass 


0.75 


VQCUT 


Quark virtuality cutoff (added to 


0.48 




quark masses in parton showers) 




VGCUT 


Gluon virtuality cutoff (added to 


0.10 




effective mass in parton showers) 




VPCUT 


Photon virtuality cutoff 


0.40 


CLMAX 


Maximum cluster mass parameter 


3.35 


CLPOW 


Power in maximum cluster mass 


2.00 


PSPLT(l) 


Split cluster spectrum parameter 


1.00 


PSPLT(2) 


1: light cluster, 2 heavy 6-cluster 


PSPLT(l) 


QDIQK 


Maximum scale for gluon— >^diquarks 


0.00 


PDIQK 


Gluon^diquarks rate parameter 


5.00 


QSPAC 


Cutoff for spacelike evolution 


2.50 


PTRMS 


Intrinsic Pt in incoming hadrons 


0.00 



Table 18: Adjustable parameters. 
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QCDLAM can be identified at high momentum fractions {x or z) with the funda- 
mental 5-flavour QCD scale A^. However, this relation does not necessarily 
hold in other regions of phase space, since higher order corrections are not 
treated precisely enough to remove renormalisation scheme ambiguities [13]. 

RMASS(l, 2, 3, 13) are effective light quark and gluon masses used in the hadron- 
ization phase of the program. They can be set to zero provided the parton 
shower cutoffs VQCUT and VGCUT are large enough to prevent divergences (see 
below) . 

For cluster hadronization, it must be possible to split gluons into qq, i.e. 
RMASS(13) must be at least twice the lightest quark mass. Similarly it may 
be impossible for heavy-flavoured clusters to decay if RMASS(4, 5) are too low. 

VQCUT and VGCUT are needed if the quark and gluon effective masses become 
small. The condition to avoid divergences in parton showers is 

1 1 1 

< 



Qi Qj QCDL3 

for either i or j or both gluons, where Qi = RMASS(i) + VQCUT for quarks, 
RMASS(13) + VGCUT for gluons, and QCDL3 is the three-flavour QCD scale used 
internally by HERWIG. QCDL3 is obtained by matching at the b- and c-quark 
mass scales from the internal five-flavour scale 

/151-97r2\ , ^ 
QCDL5 = QCDLAM X exp ( — j /V2 = 1.109 x QCDLAM . 

Note that, in the notation of ref. [13] and section pl2| , QCDL5 = Aphys/v^ for 
five flavours. 

VPCUT is the analogous quantity for photon emission. It now defaults to 
0.4 GeV. Previous versions defaulted to a/s, switching off such emission. Re- 
sults after experimental cuts are insensitive to its exact value in the range 0.1 
to 1.0 GeV. 

CLMAX and CLPOW determine the maximum allowed mass of a cluster made from 
quarks i and j as follows 

^CLPQw ^ CLMAX'^"'^™ + (RMASS(i) + RMASS(j))™ • 

Since the cluster mass spectrum falls rapidly at high mass, results become 
insensitive to CLMAX and CLPOW at large values of CLMAX. Smaller values of 
CLPOW will increase the yield of heavier clusters (and hence of baryons) for heavy 
quarks, without affecting light quarks much. For example, the default value 
gives no 6-baryons whereas CLPOW = 1.0 makes 6-baryons/6-hadrons about 1/4. 
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• PSPLT determines the mass distribution in the cluster sphtting Cdi €£2 + Ci^ 
when Cii is above the maximum allowed mass. The masses of Ci2 and Ci^ are 
generated uniformly in M^^^^"^. As long as the number of split clusters is small, 
dependence on PSPLT is weak. 

• QDIQK greater than twice the lightest diquark mass enables non-perturbative 
gluon splitting into diquarks as well as quarks. The probability of this is 
PDIQK X dQ/Q for scales Q below QDIQK. The diquark masses are taken to 
be the sum of constituent quark masses. Thus the default value QDIQK = 
suppresses gluon diquark splitting. 

• QSPAC is the scale below which the structure functions of incoming hadrons 
are frozen and non-valence constituent partons are forced to evolve to valence 
partons, if ISPAC = 0. For ISPAC = 2, structure functions are frozen at scale 
QSPAC, but evolution continues down to the infrared cutoff. 

• PTRMS is the width of the (gaussian) intrinsic transverse momentum distribution 
of valence partons in incoming hadrons at scale QSPAC. 

In practice, the parameters that have been found most effective in fitting data are 
QCDLAM, the gluon effective mass RMASS(13), and the cluster mass parameter CLMAX. 
Note that QSPAC, PTRMS and ENSGF do not affect lepton-lepton collisions. 

The default parameter values are based on those that were found to give good 
agreement when comparing earlier versions with event shape distributions at LEP. 
However, the substantial changes in this version mean that a re-tuning of parameters 
would be very worthwhile. 

Up-to-date details of HERWIG parameter tunes can be found via the official web 
page cited in section ||. 

6. Control switches, constants and options 

A number of quantities can be reset to control the program and various options: 



Name 


Description 


Default 


NEVHEP 


Current number of events 





NHEP 


Current number of entries in /HEPEVT/ 





IPRINT 


Information to include in print out 


1 


MAXPR 


Number of events to print out 


1 


PRVTX 


Include vertex information in print out 


. TRUE . 


NPRFMT 


Controls number of sig. figs, in print out 


1 


PRNDEC 


Use decimal/hexadecimal in print out 


. TRUE . 


PRNDEF 


Produce ASCII (stout) version of print out 


. TRUE . 


PRNTEX 


Produce ETgX version of print out 


.FALSE. 



Table |19| : Control switches, constants and options. (Continues) 
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Name 


Description 


Default 


PRNWEB 


Produce html version of print out 


.FALSE. 


MAXER 


Maximum number of errors to tolerate 


10 


LWEVT 


Unit for writing output events 





LRSUD 


Unit for reading Sudakov table 





LWSUD 


Unit for writing Sudakov table 


77 


SUDORD 


as order in Sudakov table 


1 


INTER 


Order of interpolation in Sudakov tables 


3 


NRN(l) 


Random number seed 1 


17673 


NRN(2) 


Random number seed 2 


63565 


WGTMAX 


Max. weight (0 to search for it) 


0.0 


NOWGT 


Generate unweighted events with EVWGT=AVWGT 


.TRUE. 


AVWGT 


Mean event weight 


1.0 


EFFMIN 


Min. acceptable Monte Carlo efficiency 


0.001 


NEGWTS 


Whether or not to allow negative weight events 


.FALSE. 


AZSOFT 


Include soft gluon azimuthal correlations 


. TRUE . 


AZSPIN 


Include gluon spin azimuthal correlations 


. TRUE . 


HARDME 


Use hard matrix-element corrections 


. TRUE . 


SOFTME 


Use soft matrix-element corrections 


. TRUE . 


GCUTME 


Gluon energy cut in top M.E. correction 


2.0 


NCOLO 


Number of colours 


3 


NFLAV 


Number of (producible) flavours 


6 


MODPDF(I) 


PDFLIB parton set and author group for beam 


-1 


AUTPDF(I) 


I (=1,2) (if MODPDF < do not use PDFLIB) 


'MRS' 


NSTRU 


Input parton set (1,2 = Duke-Owens sets 1,2; 






3,4 = EHLQ sets 1,2; 5 = Owens set 1.1, 


8 




6,7,8 = MRST, see table D 




PRSOF 


Probability of soft underlying event 


1.0 


ENSOF 


Multiplicity enhancement for SUE: 






n = (upp) (ENSOFa/s) 


1.0 


PMBNl 


Mean multiplicity in SUE/Min. bias event 


+9.110 


PMBN2 


{npp){y/s) = PMBNls™™2 + PMBN3 


+0.115 


PMBN3 




-9.500 


PMBKl 


Negative binomial param. 






fc-i = PMBKl logg(s) + PMBK2 


+0.029 


PMBK2 




-0.014 


PMBMl 


Soft cluster mass spectrum: 






{M-M1-M2- PMBMl)e-™^"'^ 


0.2 


PMBM2 




2.0 



Table |19^ Control switches, constants and options. (Continues) 
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Name 


Description 


Default 


PMBPl 


Soft cluster Pt spectrum: ^^6"^^''^^^+^'^', 






d, u quarks 


5.2 


PMBP2 


s, c quarks 


3.0 


PMBP3 


diquarks 


5.2 


lOPREM 


Options for treatment of remnant clusters 


1 


BTCLM 


Mass parameter in remnant fragmentation 


1.0 


VMIN2 


Min. parton virtuality^ in distance calcs. 


0.1 


CLRECO 


Include colour rearrangement 


.FALSE. 


PRECO 


Probability for rearrangement 


1/9 


EXAG 


Lifetime scaling for weak bosons 


1.0 


ETAMIX 


77/77' mixing angle in degrees 


-23 


PHIMIX 


(f)/uj mixing angle in degrees 


+36 


HIMIX 


/7i(1380)//7i(1170) mixing angle in degrees 


tan-i(l/v^) 


FOMIX 


— //o(1370) mixing angle in degrees 


tan-i(l/V2) 


FIMIX 


/i ( 1420) //i( 1285) mixing angle in degrees 


tan-i(l/V2) 


F2MIX 


/g/ /2 mixing angle in degrees 


+26 


ET2MIX 


7/2 (1645) 77/2(1 870) mixing angle in degrees 


tan-i(l/V2) 


OMHMIX 


— /c(j(1600) mixing angle in degrees 


tan-i(l/v^) 


PH3MIX 


(p-i/uj^ mixing angle in degrees 


+28 


BILIM 


B cluster 1 hadron parameter 


0.0 


CLDIR(I) 


Decay orientation of perturbative clusters, 


1,1 




0: isotropic, 1: along quark direction 




CLSMR(I) 


Width of gaussian angle smearing. 


0.0,0.0 




(1=1: light cluster, 1=2: heavy 6-cluster) 




PWT(I) 


A priori weights for //-pairs in cluster decay. 


1.0 




1=1-6: f = d,u, s, c, b, t 1=7: / = qq' 




REPWT(L, J,N) 


A priori weight for n^'^^^^^Lj mesons 


1.0 


SNGWT 


A priori weight for singlet baryons 


1.0 


DECWT 


A priori weight for decuplet baryons 


1.0 


PLTCUT 


Minimum lifetime for particle to be set stable 


1.0 X 10"^ 


VTOCDK(I) 


Veto decay of clusters to hadron I 


.FALSE. 


VTOCDK(I) 


Veto decay of resonances to hadron I 


. FALSE . 




1=290-293, /o(980),ao(980) 


. TRUE . 


PIPSMR 


Smear the primary vertex 


.FALSE. 


VIPWID(l) 


X width (mm) 


0.25 


VIPWID(2) 


y width (mm) 


0.015 


VIPWID(3) 


z width (mm) 


1.8 


MAXDKL 


Veto decays outside given volume 


.FALSE. 
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Name 



Description 



Default 



lOPDKL 
DXRCYL 
DXZMAX 
DXRSPH 



Option for volume: l=cylinder, 2=spliere 
Radius for cylindrical option (mm) 
Length for cylindrical option (mm) 
Radius for spherical option (mm) 



1 

20 

500 
100 



BDECAY 

MIXING 
XMIX(l) 
XMIX(2) 
YMIX(l) 
YMIX(2) 



Controls which B Decay package is used. 

Allowed values are: 'HERW, 'EURO' or 'CLEG^ 

Include neutral B meson mixing 

AM/r for 5° 

AM/r for S° 

Ar/2r for 5° 

Ar/2r for 5° 



'HERW 

. TRUE . 
10.0 
0.7 
0.2 
0.0 



RMASS(198) 
RMASS(199) 

GAMW 
RMASS(200) 

GAMZ 
WZRFR 
MODBOS(I) 



mass 
W~ mass 

width 
mass 
Z° width 

Use W/Z rest frame for decay parton showers 
Force decay modes for weak bosons, 
see sect. 



80.42 
RMASS(198) 

2.12 
91.188 

2.495 
. TRUE . 





RMASS(201) 
lOPHIG 
GAMMAX 

ENHANC(I) 



SM Higgs mass 

Options for large Higgs mass distribution 
Limit on range of Higgs mass distribution 
Enhancement factor for SM Higgs decay mode I 



115 
3 

10. 
1.0 



RMASS(209) 
RMASS(215) 



Hypothetical 4th generation 'bottom' quark mass 
corresponding antiquark mass 



200. 
RMASS(209) 



ALPHEM 
SWEIN 

QFCHCI) 
AFCHCI, J) 
VFCHCI, J) 

ZPRIME 
RMASS(202) 

GAMZP 
VCKMCI, J) 

SCABI 



1=1-6: d,..,t 
1=10-16: e, .., Ur 
j=l: Z, J=2: Z' 



Thompson limit value of aem(O) 
Value of sin^ 9w 
Fermion electric charge 
Fermion weak axial charge 
Fermion weak vector charge 
Include a Z' in 7*/Z° processes 
Mass of the Z' 
Width of the Z' 

Cabibbo-Kobayashi-Maskawa matrix elements: 
V^i,K=l-3: u,c,t L=l-3: d,s,b 
Value of sin^ 6'cabibbo 



0.0072993 
0.2319 



see sect. 4.2.2 



. FALSE 
500. 
5.0 

0.9512 0.0488 
0.0488 0.9492 0.002 
0.002 0.998 , 

0.0488 



EPOLN(l) Electron and positron beam 0.0 

EP0LN(2) polarizations in DIS and e+e~ 0.0 

EP0LN(3) annihilation. First two cmpts are 0.0 
Table |19^ Control switches, constants and options. (Continues) 
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Name 


Description 


Default 


PPOLN(l) 


transverse and only used in e"'"e", 


0.0 


PP0LN(2) 


3rd cmpt is longitudinal, and is 


0.0 


PP0LN(3) 


+/ — 1 for fully rli/lh polarized 


0.0 


QLIM 


Upper limit on hard process scale 


10^ 


THMAX 


Max. value of thrust in IPR0C=110-116 


0.9 


Y4JT 


Min. jet separation in IPR0C=600-656 


0.01 


DURHAM 


Use Durham/ JADE algorithm in IPR0C=600-656 
Colour interferences in IPR0C=600-656: 


. TRUE . 


I0P4JT(1) 


qqgg 0: neglect, 1: extreme 3142, 2: extreme 4123 





I0P4JT(2) 


qqqq 0: neglect, 1: extreme 4123, 2: extreme 2143 





BGSHAT 


Boson-gluon fusion scale (see below) 


. TRUE . 


BREIT 


Use Breit frame for DIS kinematics 


. TRUE . 


USECMF 


Use hadron-hadron cm.. 


. TRUE . 


NOSPAC 


Switch off spacelike showers 


.FALSE. 


ISP AC 


Changes meaning of QSPAC 
(see the earlier notes on QSPAC) 





TMNISR 


Min. value oi s/S for photon ISR 


10-^ 


ZMXISR 


Max. momentum fraction for photon ISR 


1 - 10-^ 


ASFIXD 


Values of fixed and uj = 12 logg(2)Q;5/7r 


0.25 


OMEGA 


for Mueller- Tang cross section 


0.3 


lAPHIG 


Approx. used in Higgs+jet production 
IPR0C=2300-2312 


1 


PHOMAS 


Damp structure functions for off mass-shell 
photons (0 for no damping) 


0.0 


PRESPL 


Preserve longitudinal momentum of hard cm. 


. TRUE . 


PTMIN 


Min. pt in hadronic jet production 


10.0 


PTMAX 


Max. pt in hadronic jet production 


108 


PTPOW 


^/Pt^°^ for jet sampling 


4.0 


YJMIN 


Min. jet rapidity 


-8.0 


YJMAX 


Max. jet rapidity 


+8.0 


EMMIN 


Min. dilepton mass in Drell-Yan 


10.0 


EMMAX 


Max. dilepton mass in Drell-Yan 


10^ 


EMPOW 


l/mF^°^ for Drell-Yan sampling 


4.0 


Q2MIN 


Min. in deep inelastic scattering 





Q2MAX 


Max. in deep inelastic scattering 


10^0 


Q2P0W 


1/Q2Q2P0W for j3jg sampling 


2.5 


YBMIN 


Min. and max. Bjorken-y in DIS 


0.0 


YBMAX 




1.0 
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Name 


Description 


Default 


WHMIN 


Mill, hadronic mass in 
7-induced processes (inc. DIS) 


0.0 


ZJMAX 


Max. z in J/ip production 


0.9 


Q2WWMN 
Q2WWMX 
YWWMIN 
YWWMAX 


Min. and max. in 

Equivalent Photon Approximation 

Min. and max. photon hght-cone fraction 

in Equiv. Photon Approx. 


0.0 
4.0 
0.0 
1.0 


CSPEED 
GEV2NB 


Speed of hght in vacuum (mm/s) 
Value of [hc/eY 


2.99792 X 10" 
389 379 


IBSH 
IBRN(l) 
IBRN(2) 


Number of shots for initial max. weight search 
1st random number seed for max. weight search 
2nd random number seed for max. weight search 


10 000 
1246579 
8447766 


NQEV 

ZBINM 
NZBIN 


Number of entries in Sudakov FF 
look-up table 

Max. bin size for z in spacelike branching 
Max. number of z bins in spacelike branching 


1024 
0.05 

100 


NBTRY 
NCTRY 
NETRY 
NSTRY 


Max. number of attempts to branch a parton 
Max. number of attempts to decay a cluster 
Max. number of attempts to generate a mass 
Max. number of attempts at soft subprocess 


200 
200 
200 
200 


ACCUR 


Precision for soft gaussian integration 


10-6 


RPARTY 


R-parity conservation in SUSY 


. TRUE . 


SUSYIN 
LRSUSY 


Check to see if SUSY data are already loaded 
Unit for reading SUSY data (if needed) 


.FALSE. 

66 


SYSPIN 
THREEB 
FOURB 


Spin correlations in decays 
SUSY three body decays 
SUSY four body decays 


. TRUE . 
. TRUE . 
.FALSE. 


TAUDEC 


Tau decay package (HERWIG or TAUOLA) 


HERWIG 


LHSOFT 
LHGLSF 


Generation of soft event for Les Houches interface 
Self-connected gluons for Les Houches interface 


. TRUE . 
.FALSE. 


OPTM 
lOPSTP 
lOPSH 


Optimisation of phase space 

Number of steps for phase space optimisation 

Number of weights for phase space optimisation 


. FALSE . 
10 

1000 



Table 19: Control switches, constants and options. 



Printout options are listed in table [2^. 

The contents of /HEPEVT/ can by printed by calling HWUEPR, those of /HWPART/ 
(the last parton shower) by calling HWUBPR. The logical variable PRNDEC (default 
.TRUE, unless NMXHEP > 9999) causes track numbers in event listings to be printed 
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IPRINT = Print program title only 

1 Print selected input parameters 

2 1 + table of particle codes and properties 

3 2 + tables of Sudakov form factors 



Table 20: Printout options. 

in decimal, or hexadecimal if false. The latter is necessary for very large events such 
as those generated by the HERBVI package (see above). 

The maximum number of errors MAXER refers to errors from which the program 
cannot recover without killing an event and starting a new one. Such errors are not 
necessarily a cause for grave concern because the phase space for backward evolution 
of initial-state showers is complicated and the program may occasionally step outside 
it (in which case the event weight should be zero anyway). When generating large 
numbers of events, it is advisable to increase MAXER in proportion, e.g. to MAXEV/100. 



See section |83 on form factors for details of LRSUD, LWSUD and SUDORD. 

The parameter EFFMIN sets the minimum allowed efficiency for the generation 
of unweighted events. A warning is printed once in every 10/EFFMIN weights if 
the efficiency is below 10 x EFFMIN, and running is stopped if the efficiency is below 



EFFMIN. See sect. |3] for details of NEGWTS. 

Variables HARDME and SOFTME invoke hard and soft matrix-element corrections 
respectively, as described in subsection |3.2.3| . 



If BGSHAT is . TRUE . , the scale used for heavy quark production via boson-gluon 
fusion in lepton-hadron collisions will be the hard subprocess cm. energy s. If it is 



. FALSE . , the scale used will be 



2 stu 



except in the case of J/ip + g production, where u is used. 

If BREIT is true, the kinematic reconstruction of deep inelastic events takes place 
in the Breit frame (i.e. the frame where the exchanged boson is purely spacelike, 
and coUinear with the incoming hadron). In fact the reconstruction procedure is 
invariant under longitudinal boosts, so any frame in which the boson and hadron are 
collinear would be equivalent, and it is only the transverse part of the boost that has 
an effect. The BREIT frame option becomes very inaccurate for very small It IS 
therefore only used if > 10~^ (the lab and Breit frames are anyway equivalent for 
such small Q"^). If BREIT is false, reconstruction takes place in the lab frame. 

If USECMF is true, the entire event record is boosted to the hadron-hadron cm. 
frame before event processing, and boosted back afterwards. This means that fixed- 
target simulation can be done in the lab frame, i.e. with PBEAM2 = 0. For hadronic 
processes with lepton beams, this boosting is always done, regardless of the value of 
USECMF. 
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In version 6.5, a new logical input variable, PRESPL [.TRUE.], has been intro- 
duced to control whether the longitudinal momentum (PRESPL = . TRUE . ) , or rapid- 
ity (PRESPL = .FALSE.), of the hard process centre-of-mass is preserved in hadron 
collisions after initial-state parton showering. At present the only function of this 
variable is to allow users to study the effects of momentum reshuffling, which is nec- 
essary after showering to compensate for jet masses. In future, it is anticipated that 
setting PRESPL= . FALSE . will simplify the treatment of other processes in MCONLO 
(see sect. Ot). 



The quantities from PTMIN to ZJMAX control the region of phase space in which 



events are generated and importance sampling inside those regions. See section |8.3.2 
on event weights for further details on these quantities and the use of WGTMAX and 
NOWGT. 

If hadronic processes with lepton beams are requested, the photon emission ver- 
tex includes the full transverse-momentum-dependent kinematics (the Equivalent 
Photon Approximation). The variables Q2WWMN and Q2WWMX set the minimum and 
maximum virtualities generated respectively. For normal simulation, Q2WWMN should 
be zero, and Q2WWMX should be the largest through which the lepton can be scat- 
tered without being detected. The variables YWWMIN and YWWMAX control the range 
of lightcone momentum fraction generated. 

In addition there are options to give different weights to the various flavours of 
quarks and diquarks, and to resonances of different spins. So far, these options have 
not been used. See the comments in the initialisation routine HWIGIN for details. 

6.1 Negative weights option 

In a number of new applications such as MCONLO (see sect. p.6|) , parton conflg- 
urations with negative weight are used to produce more accurate predictions, and 
therefore the possibility of negative event weights has to be considered. In general, 
a Monte Carlo program generates weights {wi} such that the estimated cross 
section is 

a = —^w^ = w. (6.1) 
The corresponding error depends on the width of the weight distribution: 

(6.2) 



If only positive weights are generated, and there exists a maximum weight Wmax, 
then unweighted events can be generated by 'hit and miss': w^' = or 1 with proba- 
bility P{w- = 1) = Wi/wmax- Then 



(6.3) 



a V N^W V A^e^max ViVe 
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where = N^w/ Wmax is the number of events generated. The time needed (espe- 
cially for detector simulation) depends mainly on the number of events. Hence the 
inefficiency of 'hit and miss' is not necessarily a disaster. This is the usual approach 
adopted in Monte Carlo event generators. 

Negative weights can be generated by subtraction procedures for matrix element 
corrections. These are not a problem of principle but prevent naive 'hit and miss'. 
To generalize 'hit and miss', one can generate unweighted events {w[ = 1) and 
'antievents' {w[ = —1) with 

sign(w-) = sign(wi) , (6.4) 
P{w[ = ±1) = \wi\/\wU^ . (6.5) 

Then 



Sa \w\\w\mF,^ — w 1 \w 



a V N^w JWl w 



(6.6) 



where = A^u,|u^|/|u^|max is now the total number of events+ antievents generated. 
Again, the time needed is almost proportional to N^., so this is tolerable as long as 
|w| ~ wJ. The cross section after any cuts that may be applied is 

a, = M(iv^_iV_) (6.7) 

where iV_|_ events and antievents pass the cuts. 

To allow for the possibility of negative weights, a new logical parameter NEGWTS 
has been introduced. The default (NEGWTS= . FALSE . ) is as before: negative weights 
are forbidden. If one is detected, a non-fatal warning is issued and the event weight 
is set to zero. 

If NEGWTS= . TRUE . , negative weights are allowed. Statistics are computed and 
printed accordingly. If unweighted events are requested (NOWGT= . TRUE . ) , the initial 
search stores the maximum and mean absolute weights, |w|max and \w\. Events and 



antievents are selected according to Eqs. (|6.4|J6.5| ) and EVWGT is reset to |t(7|sign(m 



so that the numerator in Eq. (|6.7|) is the sum of EVWGTs for contributing (anti)events. 



7. Particle data 

From HERWIG version 5.9 onwards, new 8-character particle names have been intro- 
duced and the revised 7 digit PDG numbering scheme, as advocated in the LEP2 
report [26] , has been adopted. All hadron and lepton masses are given to five signif- 
icant figures whenever possible. 

Unstable hadrons from clusters produced in both the hard and soft components 
of the event decay according to simplified decay schemes, which can be tabulated 
by specifying the print option IPRINT = 2. Decays modes are 'invented' where 



71 



IDHW 


RNAME 


IDPDG 


IDHW 


RNAME 


IDPDG 


395 


0MEGA_3 


227 


396 


PHI_3 


337 


397 


ETA_2(L) 


10225 


398 


ETA_2(H) 


10335 


399 


OMEGA (H) 


30223 








Table 21: New isoscalar states. 


IDHW 


RNAME 


IDPDG 


IDHW 


RNAME 


IDPDG 


57 


FH_1 


20333 








293 


FOPO 


9010221 


294 


FH_00 


10221 


62 


A_0(H)0 


10111 


290 


A_00 


9000111 


63 


A_0(H) + 


10211 


291 


A_0+ 


9000211 


64 


A_0(H)- 


-10211 


292 


A_0- 


-9000211 



Table 22: Re-identified/replaced states. 

necessary to make the branching ratios add up to 100%. Phase space distributions 
are assumed except where stated otherwise. See section ^]3|for the treatment of heavy 
quark decays. After a t, partonic b or quarkonium decay, secondary parton showers 
are produced by outgoing partons as discussed in ref. [11]; these are hadronized in 
the same way as primary jets. 

There have been a number of additions/changes to the defauh hadrons included 
via HWUDAT. Here the identification of hadrons follows the PDG [38, table 13.2], 
numbered according to their section 30. 

All S and P wave mesons are present including the ^Pq and ^Pi states and many 
new, excited B**, and quarkonium states. Also all D wave kaons and some 'light' 
1=3 states [7r2,p(1700) and ps]. All the baryons (singlet/octet/decuplet) containing 
up to one heavy (c, b) quark are included. 

New isoscalars states have been added to try to complete the l^-Da, l^D2 and 
l^Di muhiplets (table 

Also the states in table ^ have been rede- 
fined. The /i(1420) state completely replaces the 
/i(1520) in the I^Pq muhiplet, taking over 57. 
The /o(1370) (294) replaces the /o(980) (293) in 
the I^Pq multiplet; the latter is retained as it ap- 
pears in the decays of several other states. The 
new ao(1450) states (62-64) replace the three old 
ao(980) states (290-292) in the I^Pq muhiplet; the 
latter are kept as they appear in /i(1285) decays. 

By default production of the /o(980) and ao(980) states in cluster decays is 
vetoed. 

The mixing angles (in degrees) of all the light, 1=0 mesons can now be set using 
table 53. 



ETAMIX 


v/v', 


PHIMIX 


Uj/(f), 


HIMIX 


/ii(1170)//ii(1380). 


FOMIX 


/o(1300)//o(980). 


FIMIX 


/i(1285)//i(1510). 


F2MIX 


f2/fl 



Table 23: Mixing angles. 
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There were previously some inconsistencies and ambiguities in our conventions 



for the mixing of flavour 'octet' and 'singlet' mesons. They are now as in table 24 



Muitipiet 


Uctet 


Singlet 


Mixing 


Angle 


i Do 


V 


vj 


T7TAMTY 
til AfllA 


oo 

— — zo. 


1 3q 





ijj 


rniMiA 


= + O0. 


1 IP 






141 MT Y — 

nlrllA — 


illMuJ-i El 




missing 


/o(1370) 


FOMIX= 


ANGLE 


l3Pi 


/i(1420) 


/i(1285) 


F1MIX= 


ANGLE 




fl 


h 


F2MIX= 


=+26. 


1^D2 


r/2(1645) 


r?2(1870) 


ET2MIX= 


= ANGLE 


l^Di 


missing 


cj(1600) 


OMHMIX= 


= ANGLE 


1^D3 


03 




PH3MIX 


=+28. 



Table 24: Mixed meson states. 

After mixing, the quark content of the physical states is given in terms of the 
mixing angle, ^, by table |25| where tan^o = \/2. 
Hence, using the default value of 



State 


{dd + 


uu)/V2 


ss 


Octet 


cos(/ 




- sm{e + ^o) 


Singlet 


sin(^ 


^ + 0o) 


cos(6' + 6o) 



Table 25: Quark content of mixed states. 



ANGLE = arctan(l/v/2) = +35.3 for 6 
gives ideal mixing, that is, the 'octet' 
state= ss and the 'singlet'= {dd + 
uu)/\/2. This choice is important to 
avoid large isospin violations in the I^Pq 
and l^Di multiplets in which the octet member is unknown. 

Since version 6 contains a large number of supersymmetry processes many new 
hypothetical particles have been added — see section |3.5| . These new states do not 
interfere with the user's ability to add new particles as described below. 

The resonance decay tables supplied in the program have also been largely re- 
vised. Measured/expected modes with branching fraction at or above 1 per mille are 
given, including 4 and 5 body decays. To print the new tables call HWUDPR. 

The layout of HWUDAT has been altered to make it easier to identify and modify 
particle properties. Three new arrays have been introduced RLTIM, RSPIN and IFLAV. 
These are: the particle's lifetime (s), spin, and a code which specifies the flavour 
content of each hadron — used (in HWURES) to create sets of iso-flavour hadrons for 
cluster decay. Using the standard numbering of quark flavours the convention is: 

• Mesons: n^ng, e.g. vr"^: 21, vr^: 12; 

• Baryons: ±nging2^^g3, e.g. S°: 332, S°: —332 etc. (< for antibaryons; digits 
in decreasing order); 

• Light, neutral mesons are identified as: 11 if isovector (vr^, p°, . . .), 33 if isoscalar 
(r/,r/', ...). 
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Some parts of the program have been automated so that it is possible for the 
user to add new particles by specifying their properties via the arrays in /HWPROP/ 
and /HWUNAM/ and increasing NRES appropriately: this should be done before a call 
to HWUINC. 

As an example, the following lines add an isoscalar, spin-2 state 'STAN' and a 
(very light) stable toponium state 'BEER' with the decay mode: STAN BEER + 
BEER + BEER. 

NRES=NRES+1 
RNAME (NRES) =' STAN 
IDPDG (NRES) =666 
IFLAV(NRES)=11 
ICHRG(NRES)=0. 
RMASS (NRES) =0.5 
RLTIM(NRES)=1.000D-10 
RSPIN (NRES) =2.0 
NRES=NRES+1 
RNAME (NRES) =' BEER 
IDPDG (NRES) =66 
IFLAV(NRES)=66 
ICHRG(NRES)=0. 
RMASS (NRES) =0.1 
RLTIM(NRES)=1.000D+30 
RSPIN (NRES) =0.0 

CALL HWMODK (666 , 1 . DO , , 66 , 66 , 66 , , 0) 

Using the logical arrays VTOCDK and VTORDK the production of specified particles 
can be stopped in both cluster decays and via the decay of other unstable resonances. 

A priori weights for the relative production rates in cluster decays of mesons and 
baryons differing only via their S and L quantum numbers can be supplied using 
SNGWT and DECWT for singlet (i.e. A-like) and decuplet baryons and REPWT for mesons. 
The old VECWT now corresponds to REPWT(0,1,0) and TENWT to REPWT(0,2,0). 

The arrays FBTM, FTOP and FHVY which stored the branching fractions of the 
bottom, top and heavier quarks' 'partonic' decays are now no longer used. Such 
decays are specified in the same way as all other decay modes: this permits different 
decays to be given to individual heavy hadrons. Partonic decays of charm hadrons 
and quarkonium states are also now supported. As already mentioned, the products' 
order in a partonic decay mode is significant: see discussion in section |3.3| . 

The structure of the program has been altered so that the secondary hard sub- 
process and subsequent fragmentation associated with each partonic heavy hadron 
decay appears separately. Thus pre-hadronization top quark decays are treated in- 
dividually, as are any subsequent bottom hadron partonic decays. 
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Additionally decays of heavy hadrons to exclusive non-partonic final states arc 
supported. No check against double counting from partonic modes is included. How- 
ever this is not expected to be a major problem for the semi-leptonic and 2-body 
hadronic modes supplied. 

B decays can also be performed by the EURODEC or CLEO Monte Carlo pack- 
ages. The new variable BDECAY controls which package is used: 'HERW for HERWIG; 
'EURO' for EURODEC; 'CLEO' for CLEO. The EURODEC package can be obtained 
from the CERN library. The CLEO package is available by kind permission of the 
CLEO collaboration. 

An array NME has been introduced to enable a possible matrix element to be 
specified for each decay mode. 

• NME = 0: Isotropic decay. 

• NME = 100: Free particle {V -A)*{V -A), (po • P2)(pi • Ps)- 

• NME = 101: Bound quark (V - A) * (V - A), (po • P2)(Pi • ba ~ ^sPo]), = 
ttiq/Mo = spectator quark momentum fraction. 

• NME = 130: Ore and Powell ortho-positronium matrix element for: onium— > 
99 + 9h- 

• NME = 200: Free-particle i — > 6 quark decay through a scalar-fermion-fermion 
current. 

• NME = 300: Gaugino and gluino three-body decays. This also implements 
the angular ordering procedure in the gluino decays. 

The list of matrix elements currently supported is modest; users are urged to 
contact an author to have others implemented. It should be noted however that 
a number of additional matrix elements are now available via the spin correlation 
algorithm. 

A Z' has been introduced with PDG code 32, HERWIG identifier 202, default 
mass 500 GeV, width GAMZP (default 5GeV) and name 'ZOPR '. It is invoked by 
setting ZPRIME=. TRUE, (default .FALSE.). 

The decay tables can be written to/read from a file by using HWIODK, adopting 
the format advocated in the LEP2 report [26]. In addition to the PDG numbering 
of particles the HERWIG numbers or character names can be used. This permits 
easy alteration of the decay tables. In HWUINC a call is made to HWUDKS which sets 
up HERWIG internal pointers and performs some basic checks of the decay tables. 
Each decay mode must conserve charge and be kinematically allowed and not contain 
vetoed decay products. The sum of all branching ratios is set to 1 for all particles. 
Also a warning is printed if an antiparticle does not have all the charge conjugate 
decays modes of the particle. 
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HWMODK enables changes to the decay tables to be made by altering or adding 
single decay modes including on an event-by-event basis. This can be done before 
calling HWUINC, in which case when altering the branching ratio and/or matrix ele- 
ment code of an existing mode a warning is given of a duplicate second mode which 
superseeds the first. Branching ratios set below 10~^ are eliminated, whilst if one 
mode is within 10~^ of unity all other modes are removed. Note that some fore- 
thought is required if the branching ratios of two modes of the same particle are 
changed since the operation of rescaling the branching ratio sum to unity causes a 
non-commutativity in the order of the calls. 

It is possible to create particle property and event listings in any combination 
of 3 formats — standard ASCII, or html. These options are controlled by 

the logical variables PRNDEF (default .TRUE.) PRNTEX (default .FALSE.) and PRNWEB 
(default .FALSE.). The ASCII output is directed to stout (screen/log file) as in 
previous versions. When a listing of particle properties is requested (IPRINT.GE.2 
or HWUDPR is called explicitly) then the following files are produced: 

If (PRNTEX): HW_decays.tex 
If (PRNWEB): HW_decays/index.html 
/PARTOOOOOOl.html etc. 

The HW_decays .tex file is written to the working directory whilst the many ** .html 
files appear in the sub-directory HW_decays/ which must have been created previ- 
ously. Paper sizes and offsets for the output are stored at the top of the block 
data file HWUDAT: they may need modifying to suit a particular printer. When event 
hstings are requested (NEVHEP . LE . MAXPR or HWUEPR is called exphcitly) the following 
files are created in the current working directory: 

If (PRNTEX): HWEV_******* . tex 
If (PRNWEB): HWEV_* ****** .html 

where *******=0000001 etc. is the event number. 

Note that the html file automatically makes links to the index. html file of 
particle properties, assumed to be in the HW_decays sub-directory. 

A new integer variable NPRFMT (default 1) has been introduced to control how 
many significant figures are shown in each of the 3 event outputs. Basically NPRFMT=1 
gives short compact outputs whilst NPRFMT=2 gives long formats. 

Note that all the files use the package longtable . sty to format the tables. 

Also if NPRFMT=2 or PRVTX=.TRUE. then the WT^ files are designed to be printed in 
landscape mode. 
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8. Structure and output 



8.1 Main program 

The main program HWIGPR has the following form: 

PROGRAM HWIGPR 
C COMMON BLOCKS ARE INCLUDED AS FILE HERWIG65.INC 

INCLUDE 'HERWIG65.INC' 

INTEGER N 

EXTERNAL HWUDAT 
C— MAX NUMBER OF EVENTS THIS RUN 

MAXEV=100 
C BEAM PARTICLES 

PART1='P' 

PART2='P' 
C BEAM MOMENTA 

PBEAM1=7000. 

PBEAM2=PBEAM1 
C PROCESS 

IPR0C=3000 
C INITIALISE OTHER COMMON BLOCKS 

CALL HWIGIN 
C— USER CAN RESET PARAMETERS AT 
C THIS POINT, OTHERWISE DEFAULT 
C VALUES IN HWIGIN WILL BE USED. 

PRVTX=. FALSE. 

MAXER=MAXEV/100 

MAXPR=0 

PTMIN=100. 

C N.B. TO READ SUDAKOV FORM FACTOR FILE ON UNIT 77 
C INSERT THE FOLLOWING TWO LINES IN SUBSEQUENT RUNS 
C LRSUD=77 
C LWSUD=0 

C— READ IN SUSY INPUT FILE, IN THIS CASE LHC SUGRA POINT 2 

OPEN (UNIT=LRSUSY , FORM= ' FORMATTED ' , STATUS= ' OLD ' , ERR=999 , 
& FILE='sugra_pt2.1200.inO 

CALL HWISSP 

CLOSE (UNIT=LRSUSY) 
C— COMPUTE PARAMETER-DEPENDENT CONSTANTS 

CALL HWUINC 
C CALL HWUSTA TO MAKE ANY PARTICLE STABLE 

CALL HWUSTA ('PIO ') 
C USER'S INITIAL CALCULATIONS 

CALL HWABEG 
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C INITIALISE ELEMENTARY PROCESS 

CALL HWEINI 
C LOOP OVER EVENTS 

DO 100 N=1,MAXEV 
C INITIALISE EVENT 

CALL HWUINE 
C— GENERATE HARD SUBPROCESS 

CALL HWEPRO 
C GENERATE PARTON CASCADES 

CALL HWBGEN 
C DO HEAVY OBJECT DECAYS 

CALL HWDHOB 
C— DO CLUSTER FORMATION 

CALL HWCFOR 
C DO CLUSTER DECAYS 

CALL HWCDEC 
C DO UNSTABLE PARTICLE DECAYS 

CALL HWDHAD 
C— DO HEAVY FLAVOUR HADRON DECAYS 

CALL HWDHVY 
C ADD SOFT UNDERLYING EVENT IF NEEDED 

CALL HWMEVT 
C FINISH EVENT 

CALL HWUFNE 
C— USER'S EVENT ANALYSIS 

CALL HWANAL 
100 CONTINUE 
C TERMINATE ELEMENTARY PROCESS 

CALL HWEFIN 
C— USER'S TERMINAL CALCULATIONS 

CALL HWAEND 

STOP 

999 WRITE (6,*) 

WRITE (6,*) 'SUSY input file did not open correctly.' 

WRITE (6,*) 'Please check that it is in the right place.' 

WRITE (6,*) 'Examples can be obtained from the ISAWIG web page.' 

WRITE (6,*) 

END 

The declaration EXTERNAL HWUDAT is recommended to help the hnker with find- 
ing the block data on some systems. 

Various phases of the simulation can be suppressed by deleting the corresponding 
subroutine calls, or different subroutines may be substituted. For example, in non- 
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SUSY studies the call to HWISSP should be omitted, and in studies at the parton 
level everything from CALL HWCFOR to CALL HWMEVT can be omitted. 

Note that the functionality of the routine HWUINE in earlier versions has now 
been split between it and a new routine, HWUFNE. A call to the latter must be made 
between the calls to HWMEVT and HWANAL. A check is built in to prevent execution if 
this is not done. 

The analysis routine HWANAL should always begin with the line 

IF (lERROR.NE.O) RETURN 

since if an event is cancelled, each of the routines is still called in turn until reaching 
the end of the main loop. 

8.2 Form factor file 

HERWIG uses look-up tables of Sudakov form factors for the evolution of initial- and 
final-state parton showers. These can be read from an input file rather than being 
recomputed each time. The reading, writing and computing of form factor tables is 
controlled by integer parameters LRSUD and LWSUD (table pB]) . 



LRSUD 


= N > 


Read form factors for this run from unit N 


LRSUD 


= 


Compute new form factor tables for this run 


LRSUD < 


Form factor tables are already loaded 


LWSUD 


= N > 


Write form factors on unit N for future use 


LWSUD 


= 


Do not write new form factor tables 



Table 26: Form factor read/write options. 



The option LRSUD < allows the program to be initialised several times in the 
same run (e.g. to generate various event types) without recomputing or rereading 
form factors. 

Note that the Sudakov form factors depend on the parameters QCDLAM, VQCUT, 
VGCUT, NCOLO, NFLAV, RMASS(13) and RMASS(z) for i = 1,...,NFLAV. Consequently 
form factor tables must be recomputed every time any of these parameters is changed. 
These parameters are written/read with the form factor tables and checks are per- 
formed to ensure consistency. 

The parton showering algorithm uses the two-loop running coupling, with match- 
ing at each flavour threshold. However, the Sudakov table can be computed with 
either the one- loop or two- loop form, according to the variable SUDORD (= 1 or 2 re- 
spectively, default=l). If SUDORD = 1 the two-loop value is recovered using the veto 
algorithm in the shower, whereas if SUDORD = 2 no vetoes are used in the final-state 
evolution. This means that the relative weight of any shower configuration can be 
calculated in a closed form, hence that showers can be 'forced'. 
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To next-to-leading order the two possibilities should be identical, but they differ 
at beyond-NLO, so some results may change a little. The most noticeable difference 
is that the form factor table takes a factor of about five times longer to compute 
with SUDORD = 2 than with 1. 

When SUD0RD=2, no veto is needed for gluon splitting to quarks. This means that 
no vetoes are needed for final state showering, except for the previously-mentioned 
transverse momentum cut. The removal of vetoes allows preselection of the flavours 
that a jet will contain, giving a huge increase in the efficiency of rare process simu- 
lation [87]. 

8.3 Event data 

/HEPEVT/ is the LEP standard common block containing current event data (ta- 
ble ETl). 



NEVHEP 
NHEP 

ISTHEP(I) 
IDHEP(I) 
JM0HEP(1, I) 
JM0HEP(2, I) 
JDAHEP(1, I) 
JDAHEP(2, I) 
PHEP(*, I) 
VHEPf*, I) 



event number 

number of entries for this event 

status of entry I (see below) 

identity of entry I (Particle Data Group code) 

pointer to first mother of entry I (see below) 

pointer to second mother of entry I (see below) 

pointer to first daughter of entry I (see below) 

pointer to last daughter of entry I (see below) 

A/abs( 



m 



iPx,Py,Pz, E, M) of entry I: M = sign 
{x,y,z,t) of production vertex of entry I (see sect. 



Table 27: Current event data. 



All momenta are given in GeV/c in the laboratory frame, in which the input 
beam momenta are PBEAMl and PBEAM2 as specified by the user and point along the 
+z and —z directions respectively. Final state particles have ISTHEP(l) = 1. See 
section [8.3. 1| for a complete list of the special status codes used by HERWIG. 

The identity codes IDHEP are almost as recommended by the LEP Working 
Group [26], i.e. the revised PDG [38] numbers where defined, IDHEP = 91 for clus- 
ters, 94 for jets, and for objects with no PDG code. The only exception is our 
use of IDHEP = 26 for the lightest MSSM Higgs boson, to distinguish it from the 
SM Higgs boson (PDG code 25). In addition, the 'generator-specific' (IDHEP=81- 
100) codes 98 and 99 are used for remnant photons and nucleons, respectively (see 
section |3.7.2[ ) . 

HERWIG also has its own internal identity codes IDHW(I), stored in /HWEVNT/. 
The utility subroutine HWUIDT translates between HERWIG and PDG identity codes. 
See section [4.13| for further details. 
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The mother and daughter pointers are standard, except that JM0HEP(2, 1) and 
JDAHEP(2, 1) for a parton are its colour mother and colour daughter, i.e. the partons 
to which its colour and anticolour are connected, respectively. For this purpose the 
primary partons from a hard subprocess are all regarded as outgoing (see examples 
in section |4.6.2| , |4.13| and ^.5.1]) . Since a quark has no anticolour, JDAHEP(2, 1) is used 
to point to its flavour partner. Similarly for JM0HEP(2, 1) in the case of an antiquark. 

In addition to entries representing partons, particles, clusters etc, /HEPEVT/ con- 
tains purely informational entries representing the total centre-of-mass momentum, 
hard and soft subprocess momenta, etc. See section |8. 3. 1| for the corresponding status 
codes. 

Information from all stages of event processing is retained in /HEPEVT/ so the 
same particle may appear several times with different status codes. For example, 
an outgoing parton from a hard scattering (entered initially with status 113 or 114) 
will appear after processing as an on-mass-shell parton before QCD branching (sta- 
tus 123,124), an off-mass-shell entry representing the flavour and momentum of the 
outgoing jet (status 143,144), and a jet constituent (157). It might also appear again 
in other contexts, e.g. as a spectator in a heavy flavour decay (status 154,160). 

Incoming partons (entered with status 111 or 112, changed to 121 or 122 after 
branching) give rise to spacelike jets (status 141 or 142) with < 0, indicated by 
PHEP(5, IHEP) < 0, due to the loss of momentum via initial-state bremsstrahlung. 
The same applies in principle to incoming leptons, but in that case emission of at 
most a single photon is permitted from each initial-state lepton and the off-shell 
lepton is given status code 3 (see section p.3.1| ). 

Each parton jet begins with a status 141-144 jet entry giving the total flavour 
and momentum of the jet. The first mother pointer of this entry gives the location 
of the parent hard parton, while the second gives that of the subprocess centre-of- 
mass momentum. If QCD branching has occurred, this is followed by a lightlike 
CONE entry, which fixes the angular extent of the jet and its azimuthal orientation 
relative to the parton with which it interferes. The interfering parton is listed as 
the second mother of the cone. Next come the actual constituents of the jet. If no 
branching has occurred, there is no cone and the single jet constituent is the same 
as the jet. 

Since version 5.1, the event record has been modified to retain entries for all 
partons before hadronization (with status ISTHEP=2). During hadronization, the 
gluons are split into quark-antiquark, while other partons are copied to a location 
(indicated by JDAHEP(1,*)) where their momenta may be shifted slightly, to conserve 
momentum, during heavy cluster splitting. Previously the original momenta were 
shifted, so momentum appeared not to be conserved at the parton level. 

Since version 6.3, to take account of the increased energy and complexity of inter- 
actions at LHC and future colliders, the default value of the parameter NMXHEP, which 
sets the array sizes in the standard /HEPEVT/ common block, has been increased to 



81 



4000. 



8.3.1 Status codes 

A complete list of currently-used HERWIG status codes is given below (table [28| ). 
Many are used only in intermediate stages of event processing. The most important 
for users are probably 1 (final-state particle), 101-3 (initial state), 141-4 (jets), and 
199 (decayed 6-flavoured hadrons). 

For technical reasons, some HERWIG status codes ISTHEP between 153 and 165 
have changed their meanings since version 5.1. 

The event status ISTAT in common /HWEVNT/ is roughly ISTHEP — 100 where 
ISTHEP is the status of entries being processed. For completed events, ISTAT = 100. 





Description 


-1 

1 


nnal state particle 




parton before hadronization 


Q 

o 


documentation line 




cone limiting jet evolution 


iUi 


'beam' (beam 1) 


1 no 


'target' (beam 2) 


1 no 
iUo 


overall centre of mass 


110 


unprocessed hard process cm. 


111 

iii 


unprocessed beam parton 


1 1 o 

LIZ 


unprocessed target parton 


1 1 o 

113 


unproc. nrst outgomg parton 


1 14 


1 iTi TiTor' f^i" n (^r r^i i or^i no" T»fi yir^Ti 

LiiijJivJL/. CLliCi CU-L^VJllig, jJdi LCli 


115 


unprocessed spectator parton 


120-25 


as 110-15, after processing 


130 


lepton in jet (unboosted) 


131-34 


as 141-44, unboosted to cm. 


135 


spacelike parton (beam, unboosted) 


136 


spacelike parton (target, unboosted) 


137 


spectator (beam, unboosted) 


138 


spectator (target, unboosted) 


139 


parton from branching (unboosted) 


140 


parton from gluon splitting (unboosted) 


141-44 


jet from parton type 111-14 


145-50 


as 135-40 boosted, unclustered 


151 


as 159, not yet clustered 


152 


as 160, not yet clustered 


153 


spectator from beam 



Table Status codes. (Continues) 
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ISTHEP 


Description 


154 


spectator from target 


155 


heavy quark before decay 


156 


spectator before heavy decay 


157 


parton from QCD branching 


158 


parton from gluon sphtting 


159 


parton from cluster sphtting 


160 


spectator after heavy decay 


161 


beam spectator after gluon splitting 


162 


target spectator after gluon splitting 


163 


other cluster before soft process 


164 


beam cluster before soft process 


165 


target cluster before soft process 


167 


unhadronized beam cluster 


168 


unhadronized target cluster 


170 


soft process centre of mass 


171 


soft cluster (beam, unhadronized) 


172 


soft cluster (target, unhadronized) 


173 


soft cluster (other, unhadronized) 


181 


beam cluster (no soft process) 


182 


target cluster (no soft process) 


183 


hard process cluster (hadronized) 


184 


soft cluster (beam, hadronized) 


185 


soft cluster (target, hadronized) 


186 


soft cluster (other, hadronized) 


190-93 


as 195-98, before decays 


195 


direct unstable non-hadron 


196 


direct unstable hadron (1-body clus.) 


197 


direct unstable hadron (2-body clus.) 


198 


indirect unstable hadron or lepton 


199 


decayed heavy flavour hadron 


200 


neutral B meson, flavour at prod'n 



Table 28: Status codes. 



8.3.2 Event weights 

The default is to generate unweighted events (EVWGT = AVWGT). Then event distri- 
butions are generated by computing a weight proportional to the cross section and 
comparing it with a random number times the maximum weight. Set WGTMAX to the 
maximum weight, or to zero for the program to compute it. If a weight greater than 
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WGTMAX is generated during execution, a warning is printed and WGTMAX is reset. If 
this occur too often, output event distributions could be distorted. 

To generate weighted events, set NOWGT = .FALSE, in common /HWEVNT/. 

In QCD hard scattering and heavy flavour and direct photon production (IPR0C= 
1500-1800) the transverse energy distribution of weighted events (or the efficiency 
for unweighted events) can be varied using the parameters PTMIN, PTMAX and PTPOW. 

Similarly in Drell-Yan processes (IPROC = 1350 etc.) the lepton pair mass distri- 
bution is controlled by the parameters EMMIN, EMMAX and EMPOW, and in deep inelastic 
scattering the distribution depends on Q2MIN, Q2MAX and Q2P0W. 

Data on weights generated are output at the end of the run. The mean weight 
is an estimate of the cross section (in nanobarns) integrated over the region used for 
event generation. Note that the mean weight is the sum of weights divided by the 
total number of weights generated, not the total number of events. 

The maximum weight is now always printed in full precision. This is needed to 
be sure of generating the same events in repeated runs. 

In versions 6.3 and higher, the option of generating negative weights has been 
available. This is controlled by the input parameter NEGWTS [.FALSE.], see sect. ^ 
for details. 

8.4 Error conditions 

Certain combinations of input parameters may lead to problems in execution. HER- 
WIG tries to detect these and print a warning. Errors during execution are dealt with 
by HWWARN, which prints the calling subprogram and a code and takes appropriate 
action. In general, the larger the code the more serious the problem. Refer to the 
source code to find out why HWWARN was called. It is important to note the subpro- 
gram from which the call was issued: many different subprograms use the same error 
code, but each code is unique within a given subprogram. 

Events can be rerun by setting the random number seeds NRN(l) and NRN(2) to 
the values given in the error message or event dump, and MAXWGT to the maximum 
weight encountered in the run. The contents of /HEPEVT/ can by printed by calling 
HWUEPR, those of /HWPART/ (the last parton shower) by calling HWUBPR. 

Note that if WGTMAX is increased during event generation, so that this type of 
message is printed: 

HWWARN CALLED FROM SUBPROGRAM HWEPRO : CODE = 1 

EVENT 21: SEEDS = 836291635 & 1823648329 WEIGHT = 0.3893E-08 

EVENT SURVIVES. EXECUTION CONTINUES 

NEW MAXIMUM WEIGHT = . 428217360829367E-08 

then to regenerate any later events, WGTMAX must be set to the printed value, as well 
as setting NRN to the appropriate seeds. 
Examples of error messages are: 
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HWWARN CALLED FROM SUBPROGRAM HWSBRN: CODE = 101 

EVENT 31: SEEDS = 422399901 & 771980111 WEIGHT = 0.3893E-08 

EVENT KILLED. EXECUTION CONTINUES 

Spacelike (initial-state) parton branching had no phase space. This can happen due 
to cutoffs which are slightly different in the hard subprocess and the parton shower. 
Action taken: program throws away this event and starts a new one. 

HWWARN CALLED FROM SUBPROGRAM HWCHAD: CODE = 102 

EVENT 51: SEEDS = 1033784787 & 1428957533 WEIGHT = 0.3893E-08 

EVENT KILLED. EXECUTION CONTINUES 

A cluster has been formed with too low a mass to represent any hadron of the correct 
flavour, and there is no colour-connected cluster from which the necessary additional 
mass could be transferred. 

Action taken: program throws away this event and starts a new one. 

HWWARN CALLED FROM SUBPROGRAM HWUINE: CODE= 200 
EVENT SURVIVES. RUN ENDS GRACEFULLY 

CPU time limit liable to be reached before generating MAXEV events. 
Action taken: skips to terminal calculations using existing events. 

HWWARN CALLED FROM SUBPROGRAM HWBSUD: CODE= 500 
RUN CANNOT CONTINUE 

The table of Sudakov form factors read on unit LRSUD does not extend to the maxi- 
mum momentum scale QLIM specified for this run. 

Action taken: run aborted. The user must either reduce QLIM or set LRSUD to zero 
to make a bigger table (set LWSUD non-zero to write it). 

HWWARN CALLED FROM SUBPROGRAM HWBSUD: CODE= 515 
RUN CANNOT CONTINUE 

The table of Sudakov form factors read on unit LRSUD is for a different value of a 
relevant parameter (in this case the b quark mass). 

Action taken: run aborted. The user must make a new table (set LWSUD non-zero to 
write it). 

8.5 Sample output 



This is the output from the main program listed in section with no event printout 
or user analysis. 



HERWIG 6.500 16 Oct 2002 



Please reference: G. Marchesini, B.R. Webber, 
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G.Abbiendi, I.G.Knowles, M.H.Seymour & L. Stance 
Computer Physics Communications 67 (1992) 465 

and 

G.Corcella, I.G.Knowles, G.Marchesini, S.Moretti, 
K.Odagiri, P.Richardson, M.H.Seymour & B.R.Webber, 
JHEP 0101 (2001) 010 



Since SUSY processes are called, 
please also reference: S.Moretti, K.Odagiri, 
P.Richardson, M.H.Seymour & B.R.Webber, 
JHEP 0204 (2002) 028 

Reading in SUSY data from unit 66 

INPUT CONDITIONS FOR THIS RUN 



BEAM 1 (P ) MOM. = 7000.00 

BEAM 2 (P ) MOM. = 7000.00 

PROCESS CODE (IPROC) = 3000 

NUMBER OF FLAVOURS = 6 

STRUCTURE FUNCTION SET = 8 

AZIM SPIN CORRELATIONS = T 

AZIM SOFT CORRELATIONS = T 

QCD LAMBDA (GEV) = 0.1800 

DOWN QUARK MASS = 0.3200 

UP QUARK MASS = 0.3200 

STRANGE QUARK MASS = 0.5000 

CHARMED QUARK MASS = 1.5500 

BOTTOM QUARK MASS = 4.9500 

TOP QUARK MASS = 175.0000 

GLUON EFFECTIVE MASS = 0.7500 

EXTRA SHOWER CUTOFF (Q)= 0.4800 

EXTRA SHOWER CUTOFF (G)= 0.1000 

PHOTON SHOWER CUTOFF = 0.4000 

CLUSTER MASS PARAMETER = 3.3500 

SPACELIKE EVOLN CUTOFF = 2.5000 

INTRINSIC P-TRAN (RMS) = 0.0000 

DECAY SPIN CORRELATIONS= T 

SUSY THREE BODY ME = T 

SUSY FOUR BODY ME = F 
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NO EVENTS WILL BE WRITTEN TO DISK 

B_d: Delt-M/Gam =0.7000 Delt-Gam/2*Gam =0.0000 
B_s: Delt-M/Gam = 10.00 Delt-Gam/2*Gam =0.2000 

PDFLIB NOT USED FOR BEAM 1 
PDFLIB NOT USED FOR BEAM 2 



Checking consistency of particle properties 



Checking consistency of decay tables 



Line 




i 


± O liLLfs b cLIIlcj do J- lILfs 








Take 


BR 


0. 


.333 and ME code 100 


from 


second 


entry 


Line 




2 


is the same as line 


2635 






Take 


BR 


0. 


.333 and ME code 100 


from 


second 


entry 


Line 




3 


is the same as line 


2636 






Take 


BR 


0, 


.111 and ME code 100 


from 


second 


entry 


Line 




4 


is the same as line 


2637 






Take 


BR 


0. 


.111 and ME code 100 


from 


second 


entry 


Line 




5 


is the same as line 


2638 






Take 


BR 


0. 


.111 and ME code 100 


from 


second 


entry 


Line 




6 


is the same as line 


2639 






Take 


BR 


0, 


.333 and ME code 100 


from 


second 


entry 


Line 




7 


is the same as line 


2640 






Take 


BR 


0. 


.333 and ME code 100 


from 


second 


entry 


Line 




8 


is the same as line 


2641 






Take 


BR 


0. 


.111 and ME code 100 


from 


second 


entry 


Line 




9 


is the same as line 


2642 






Take 


BR 


0, 


.111 and ME code 100 


from 


second 


entry 


Line 


10 


is the same as line 


2643 






Take 


BR 


0, 


.111 and ME code 100 


from 


second 


entry 



WRITING SUDAKOV TABLE ON UNIT 77 



WRITING MATRIX ELEMENT TABLE ON UNIT 88 
CHECKING SUSY DECAY MATRIX ELEMENTS 
PARTICLE TYPE 21=PI0 SET STABLE 
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INITIAL SEARCH FOR MAX WEIGHT 



PROCESS CODE IPROC = 
RANDOM NO. SEED 1 = 
SEED 2 = 
NUMBER OF SHOTS 
NEW MAXIMUM WEIGHT = 
NEW MAXIMUM WEIGHT = 
NEW MAXIMUM WEIGHT = 



3000 
1246579 
8447766 
10000 

2 . 4359288434851706E-03 
3 . 7017893759328114E-03 
1 . 9578185739344733E-02 



HWWARN CALLED FROM SUBPROGRAM HWSMRS : CODE = 5 

EVENT 0: SEEDS = 17673 & 63565 WEIGHT = O.OOOOE+00 

EVENT SURVIVES. EXECUTION CONTINUES 

WARNING: MRST98 CALLED WITH Q OUTSIDE ALLOWED RANGE! 

Q VALUE=3 . 869E+03 , MINIMUM=1 . 118E+00, MAXIMUM=3 . 162E+03 

NO FURTHER WARNINGS WILL BE ISSUED 

NEW MAXIMUM WEIGHT = 2 . 2550328399448680E-02 
NEW MAXIMUM WEIGHT = 5 . 5633613959185035E-02 
NEW MAXIMUM WEIGHT = 6 . 9166986913034176E-02 
NEW MAXIMUM WEIGHT = 0.1197405440791255 
NEW MAXIMUM WEIGHT = 0.1488818465259009 



INITIAL SEARCH FINISHED 



OUTPUT ON ELEMENTARY PROCESS 



N.B. NEGATIVE WEIGHTS NOT ALLOWED 



NUMBER OF EVENTS 
NUMBER OF WEIGHTS 
MEAN VALUE OF WGT 
RMS SPREAD IN WGT 
ACTUAL MAX WEIGHT 
ASSUMED MAX WEIGHT 





10000 
= 4.1213E-03 
= 1.2862E-02 
= 1.3755E-01 
= 1.4888E-01 



PROCESS CODE IPROC = 3000 

CROSS SECTION (PB) = 4.121 

ERROR IN C-S (PB) = 0.1286 

EFFICIENCY PERCENT = 2.768 



SUBROUTINE TIMEL CALLED BUT NOT LINKED. 
DUMMY TIMEL WILL BE USED. DELETE DUMMY 
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AND LINK CERNLIB FOR CPU TIME REMAINING. 



OUTPUT ON ELEMENTARY PROCESS 

N.B. NEGATIVE WEIGHTS NOT ALLOWED 

NUMBER OF EVENTS = 100 
NUMBER OF WEIGHTS = 3699 

MEAN VALUE OF WGT = 3.7473E-03 

RMS SPREAD IN WGT = 1.2205E-02 

ACTUAL MAX WEIGHT = 1.3900E-01 

ASSUMED MAX WEIGHT = 1.4888E-01 

PROCESS CODE IPROC = 3000 

CROSS SECTION (PB) = 3.747 

ERROR IN C-S (PB) = 0.2007 

EFFICIENCY PERCENT = 2.517 

8.5.1 Guide to sample output 

See ref. [1] for a full discussion of the basic features of HERWIG output, including a 
listing of a sample event. Here we point out only some new features in comparison 
to version 5.1. 

The beam particles, their energies and the process code IPR0C=3000 indicate 
the SUSY process of squark and/or gluino production at LHC The call to HWISSP 



triggers a request for a SUSY particle data file, in the format specified in section 3.5 



which is read in from the default unit. In this case the file corresponds to the second 



LHC SUGRA point discussed in section p. 5.1 



First the program lists the main relevant input parameters, including and 
Bg mixing parameters. Parton distributions were not requested from the PDFLIB 
library and therefore the default MRST set [51] is used (see sect. [4.1.1| .) 

Next the program performs some basic checks on the particle data provided. 
Here it finds that the input file read by HWISSP contains top quark decay modes 
which duplicate the default modes stored in HWUDAT. The branching ratios and matrix 
element codes are accordingly updated to those in the input file. 

After an initial search for the maximum weight, the program prints its estimate 
of the relevant cross section and the expected Monte Carlo efficiency for event gen- 
eration. During the search for the maximum weight the default parton distributions 
are called with a parameter, in this case the scale, outside the allowed range and a 
warning is issued. The maximum or minimum value of the parameter is used instead, 
depending on whether the value requested is too large or too small. 

In the course of event generation a new cross section estimate is generated, which 
is printed together with the actual Monte Carlo efficiency at the end of the run. 
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8.6 Subroutine descriptions 



We give in table a list of all subroutines with their functions. Note that the third 
letter of the name usually follows a rough classification scheme. 



Name 


Description 




Main program and initialisation 


HWIGPR 


Main program 


HWIGIN 


Default initialisations 




Reading/writing/altering decay modes 


HWIGUP 


Initialisation for Les Houches interface 


HWIMDE 


Adds modes for SUSY four body decays 


HWIODK 


Inputs/outputs formatted decay tables 


HWISPC 


Initialisation of couplings for spin correlations and SUSY decays 


HWISPN 


Initialisation for spin correlations and SUSY decays 


HWISP2 


Initialisation of two body decays for spin correlations and SUSY decays 


HWISP3 


Initialisation of three body decays for spin correlations and SUSY decays 


HWISP4 


Initialisation of four body decays for spin correlations and SUSY decays 


HWIPHS 


Initialises optimized phase space 


HWISSP 


Inputs supersymmetric particle data 


HWMODK 


Modifies or adds an individual decay mode 




User-provided analysis routines 


HWABEG 


Initialises user's analysis 


HWAEND 


Terminates user's analysis 


HWANAL 


Performs user's analysis on event 




Parton branching with interfering gluons 


HWBAZF 


Computes azimuthal correlation functions 


HWBCON 


Makes colour connections between jets 


HWBDED 


Correction to the 'dead zone' in e^e~ 


HWBDIS 


Correction to the 'dead zone' in DIS 


HWBDYP 


Correction to the 'dead zone' in Drell-Yan 


HWBFIN 


Transfers external lines of jet to /HEPEVT/ 


HWBGEN 


Finds unevolved partons and generates jets 


HWBGUP 


Makes colour and flavour connections for Les Houches events 


HWBJCO 


Combines jets with correct kinematics 


HWBMAS 


Computes masses and trans, momenta in jet 


HWBRAN 


Generates a timelike parton branching 


HWBRCN 


Replaces HWBCON if R-parity is violated 


HWBRCl 


Finds colour partner in gluino decay 


HWBRC2 


Finds colour partner in jet 


HWBSPA 


Computes momenta in spacelike jet 




Table Subroutines. (Continues) 
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Name 



Description 



HWBSPN Computes spin density /decay matrices 

HWBSUl First term in quark Sudakov form factor 

HWBSU2 Second term in quark Sudakov form factor 

HWBSUD Computes (or reads) Sudakov form factors 

HWBSUG Integrand in gluon Sudakov form factor 

HWBSUL Logarithmic part of Sudakov form factor 

HWBTIM Computes momenta in timelike jet 

HWBTOP Correction to the 'dead zone' in top decay 

HWBVMC Virtual mass cutoff for parton type ID 



Cluster hadronization model 



HWCBCT Cuts a massive baryon cluster in two 

HWCBVI Clusters quarks from a p interaction 

HWCBVT Finds which p interaction partons came from 

HWCCCC Finds colour connections after gluon splitting if 

HWCCUT Cuts a massive cluster in two 

HWCDEC Decays clusters into primary hadrons 

HWCFLA Sets up flavours for HWCHAD 

HWCFOR Forms clusters 

HWCGSP Splits gluons 

HWCHAD Decays a cluster into one or two hadrons 



Particle and heavy quark decays 



HWD2ME Finds maximum weight for a two body decay for spin correlations 

HWD3ME SUSY three body decays and spin correlations master routine 

HWD3M0 Calculation of SUSY three body matrix element 

HWD3M1 Helicity amplitude for / fff via vector boson exchange 

HWD3M2 Helicity amplitude for / fff via Higgs exchange 

HWD3M3 Helicity amplitude for / fff via antisfermion exchange 

HWD3M4 Helicity amplitude for / — > /// via sfermion exchange 

HWD3M5 Helicity amplitude for / — > /// via vector boson exchange 

HWD3M6 Helicity amplitude for — > 0// via vector boson exchange 

HWD3M7 Helicity amplitude for / Gff via vector boson exchange 

HWD3M8 Helicity amplitude for / fff via scalar exchange (1st diagram) 

HWD3M9 Helicity amplitude for / fff via scalar exchange (2nd diagram) 

HWD3MA Helicity amplitude for / fff via scalar exchange (3rd diagram) 

HWD3MB Helicity amplitude for / fff via scalar exchange (1st diagram) 

HWD3MC Helicity amplitude for / — > /// via scalar exchange (2nd diagram) 

HWD3MD Helicity amplitude for / — > /// via scalar exchange (3rd diagram) 

HWD3MF Helicity amplitude for / —* fff via scalar exchange (1st diagram) 

HWD3MG Helicity amplitude for / —* fff via scalar exchange (2nd diagram) 

Table Subroutines. (Continues) 
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Name 



Description 



HWD3MH Helicity amplitude for / fff via scalar exchange (3rd diagram) 

HWD3MI Helicity amplitude for / fff via scalar exchange 

HWD4ME Master routine for SUSY four body decays 

HWD4M0 Matrix Element for (p ^ V*V* ^ ffff 

HWDBOS Finds and decays W and Z bosons 

HWDBOZ Chooses decay mode of W and Z bosons 

HWDBZ2 Copy of HWDBOZ used by hadronic WW, WZ and ZZ production 

HWDCHK Checks given decay mode is self-consistent 

HWDCLE Interface to CLEO package for B decays 

HWDEUR Interface to EURODEC package for B decays 

HWDFIV Generates a five-body decay 

HWDFOR Generates a four-body decay 

HWDHAD Generates decays of unstable hadrons 

HWDHGC Higgs 77 decay 

HWDHGF Higgs ^ W+W' decay 

HWDHIG Finds and decays Higgs bosons 

HWDHOB Finds and decays heavy objects 

HWDHOl Selects decay mode for heavy object 

HWDH02 Calculates momenta of heavy object's decay products 

HWDH03 Makes colour connects for heavy object decay 

HWDH04 Performs parton shower for heavy object decay 

HWDH05 Checks for colour disconnections in heavy object decay 

HWDH06 Perform flp colour connections in heavy object decay 

HWDHVY Finds and decays heavy flavours 

HWDHWT Subroutine for top decay via a virtual if^ 

HWDRCL Colour connections for a SUSY decay 

HWDRME Main 3-body decay matrix element subroutine 

HWDRMl fip 3-body decay matrix element subroutine 

HWDRM2 fip 3-body decay matrix element subroutine 

HWDRM3 fip 3-body decay matrix element subroutine 

HWDRM4 3-body decay matrix element subroutine 

HWDRM5 ^p 3-body decay matrix element subroutine 

HWDPWT Phase space three-body decay weight 

HWDSIN Performs decays with spin correlations 

HWDSIl Picks next particle to decay for spin correlations 

HWDSI2 Pick moment of decay products with spin correlations 

HWDSI3 Selects r polarization and passes it to the decay routine 

HWDSM2 Two body decay with spin correlations 

HWDSM3 Three body decay with spin correlations 



Table 29: Subroutines. (Continues) 
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Name 


Description 


HWDSM4 


Four body decay with spin correlations 


HWDTAU 


Interface to TAUOLA for r decays 


HWDTHR 


Generates a three-body decay 


HWDTOP 


Decides whether to decay top quark 


HWDTWO 


Generates a two-body decay 


HWDWWT 


Weak (y — A) three-body decay weight 


HWDXLM 


Tests if decay vertex hes in given volume 




Elementary subprocess generation 


HWECIR 


Interface to CIRCE 


HWEFIN 


Final calculations on elementary subprocess 


HWEGAM 


Generates incoming photon 


HWEINI 


Initialises elementary subprocess 


HWEISR 


Generates photon emission from e or /x 


HWEONE 


Sets up a 2 — i> 1 hard subprocess 


HWEPRO 


Generates elementary subprocess 


HWETWO 


Sets up a 2 — i> 2 hard subprocess 




Individual hard subprocesses 


HWH2BK 


Matrix element for bb W^H"^ 


HWH2BH 


Matrix element for production via 6g-fusion 


HWH2DD 


Function to return the D function of [68] 


HWH2F1 


Subroutine to return the F function of [69] for a fixed first momentum 


HWH2F2 


Subroutine to return the F function of [69] for a fixed second momentum 


HWH2F3 


Subroutine to return the F function of [69] for all first and second momenta 


HWH2HE 


Matrix element for Higgs associated production 


HWH2M0 


Subroutine to compute the massless matrix element for QQZ 


HWH2MQ 


Subroutine to compute the massive matrix element for QQZ 


HWH2PS 


Subroutine to perform the phase-space for Z-|-two jets 


HWH2P1 


Subroutine to select quark masses for HWH2PS 


HWH2P2 


Subroutine to select quark masses for HWH2PS 


HWH2QH 


Matrix element for qq,gg — QQ'^'^ Higgs 


HWH2SH 


Matrix element for squark pair plus Higgs production 


HWH2SS 


Subroutine to return the S function of [68] 


HWH2T1 


Function to return the Ti function of [68] 


HWH2T2 


Function to return the T2 function of [68] 


HWH2T3 


Function to return the T3 function of [68] 


HWH2T4 


Function to return the T4 function of [68] 


HWH2T5 


Function to return the T5 function of [68] 


HWH2T6 


function to return the function of [68] 


HWH2T7 


Function to return the T7 function of [68] 




Table 29: Subroutines. (Continues) 
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Name 



Description 



HWH2T8 
HWH2T9 
HWH2T0 
HWH2VH 
HWH4JT 
HWH4J1 
HWHBGF 
HWHBKI 
HWHBRN 
HWHBSG 
HWHDIS 
HWHDYP 
HWHDYQ 
HWHEGG 
HWHEGW 
HWHEGX 
HWHEPA 
HWHEPG 
HWHESG 
HWHESL 
HWHESP 
HWHESQ 
HWHEWO 
HWHEWl 
HWHEW2 
HWHEW3 
HWHEW4 
HWHEW5 
HWHEWW 
HWHGBP 
HWHGBS 
HWHGBl 
HWHGB2 
HWHGBS 
HWHGB4 
HWHGB5 
HWHGRV 
HWHGUP 
HWHHVY 



Function to return the Tg function of [68] 

Function to return the Tg function of [68] 

Function to return the Tio function of [68] 

Matrix element for qq^'^ V Higgs, V = W^, . 

Hard subprocess: e^e~ — > 4 jets 

Matrix element for e~^e~ 4 jets 

Hard subprocess: boson-gluon fusion (BGF) 

Computes kinematics for BGF 

Returns a phase-space point for BGF 

Computes cross section for BGF 

Hard subprocess: deep inelastic e//i quark 

Hard subprocess: Drell-Yan Z^/'j production 

Subroutine for QQZ 

Hard subprocess: two-photon processes in e~^e~ 
Hard subprocess: processes in e+e" 
Calculates cross section for HWHEGW 
Hard subprocess: e+e~ —>■ qq 
Hard subprocess: e+e~ — *• qqg 
Gaugino pair production in e"'"e~ collisions 
Slepton pair production in e~^e~ collisions 
Sparticle pair production in e'^e~ collisions 
Squark pair production in e'^e~ collisions 



e ' 



w+w- 
w+w- 
w+w- 
w+w- 
w+w- 
w+w- 



subroutine 
subroutine 
subroutine 
subroutine 
subroutine 
subroutine 



Hard subprocess: e+e~ W^W~ 

Main routine for gauge boson pair production in hadron-hadron 

Phase space for gauge boson pair production in hadron-hadron 

Selects gauge boson mass for HWHGBS 

Matrix element for WW in hadron-hadron 

Matrix element for ZZ in hadron-hadron 

Matrix element for WZ in hadron-hadron 

Selects t and u for HWHGBS 

Graviton resonance production 

External hard process using Les Houches interface 

Hard subprocess: heavy quark production 

Table B^: Subroutines. (Continues) 
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Name 



Description 



HWHIBG 
HWHIBK 
HWHIBQ 
HWHIGl 
HWHIGA 
HWHIGB 
HWHIGE 
HWHIGH 
HWHIGJ 
HWHIGM 
HWHIGQ 
HWHIGS 
HWHIGT 
HWHIGV 
HWHIGW 
HWHIGY 
HWHIGZ 
HWHIHH 
HWHISQ 
HWHPH2 
HWHPHO 
HWHPPB 
HWHPPE 
HWHPPH 
HWHPPM 
HWHPPT 
HWHPQS 
HWHQCD 
HWHQCP 
HWHQCM 
HWHRBB 
HWHRBS 
HWHREE 
HWHREM 
HWHREP 
HWHRES 
HWHRLL 
HWHRLS 
HWHRSP 



Hard subprocess: for bg Q Higgs, with Q = t,b 

Hard subprocess: bb — > IV^H^ 

Subroutine for production via 6g-fusion 

Matrix elements for Higgs + jet production 

Amplitudes squared for Higgs + jet 

Loop integrals for Higgs + jet 

Hard subprocess: Higgs associated production 

Hard subprocess: qq — > Higgsi + HiggS2 

QCD Higgs + jet production 

Choose any Higgs mass for production routines 

Hard subprocess: qq, gg QQ^'^ Higgs 

Hard subprocess: gg/qq-^ Higgs 

Computes gg Higgs cross section 

Hard subprocess: gg*^'^ V Higgs, V = W^, 

Hard subprocess: W^W~ jZ^Z^ Higgs 

Computes e"*"e~ Z^ ^ Z^H^-^ cross section 

Hard subprocess: £+£- ^ Z° ^ ^°^sm 

Hard subprocess: SUSY Higgs pair production 

Subroutine for squark pair plus Higgs production 

Hard subprocess: direct photon pairs 

Hard subprocess: direct photon production 

Box contribution to gg — 77 

Pointlike photon-parton (fixed flavour) 

Pointlike photon-parton (fixed pair flavour) 

Pointlike photon-parton direct light meson 

Pointlike photon-parton (all flavours) 

Pointlike photon-quark (Compton) scattering 

QCD 2^2 hard subprocesses 

Identifies QCD 2 — > 2 hard subprocess 

Hard subprocess: 77 qq, W^W~ 

fip resonant squark to SM particles 

fip resonant squark to SUSY particles 

fip SM particle production in e+e~ collisions 

Treats hard scattering remnants 

Decides which subroutine to use in e~^e~ 

fip single sparticle production in e+e~ collisions 

fip resonant slepton to SM particles 

fip resonant slepton to SUSY particles 

Decides which subroutine to use in hadron-hadron 



Table 29: Subroutines. (Continues) 
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Name 



Description 



HWHRSS 
HWHSCT 
HWHSNG 
HWHSNM 
HWHSPN 
HWHSOl 
HWHS02 
HWHS03 
HWHS04 
HWHS05 
HWHS06 
HWHS07 
HWHS08 
HWHS09 
HWHSIO 
HWHSll 
HWHS12 
HWHS13 
HWHS14 
HWHS15 
HWHS16 
HWHS17 
HWHS18 
HWHS19 
HWHS20 
HWHS21 
HWHS22 
HWHS23 
HWHS24 
HWHS25 
HWHS26 
HWHS27 
HWHS28 
HWHS29 
HWHS30 
HWHS31 
HWHS32 
HWHS33 
HWHS34 



Identifies flp process 

Process extra hard scatterings 



Co: 

Co: 

Ma: 

He 

He 

He 

He 

He 

He 

He 

He 

He 

He 

He 

He 

He 

He 

He 

He 

He 

He 

He 

He 

He 

He 

He 

He 

He 

He 

He 

He 

He 

He 

He 

He 

He 

He 



our singlet parton scattering 
our singlet parton scattering matrix element 
n routine for spin correlations in the hard process 
city amplitude for // — > // via s-channel vector boson exchange 
city amplitude for // ff via t-channel scalar exchange 
city amplitude for // ff via w-channel scalar exchange 
city amplitude for // ff via s-channel vector boson exchange 
city amplitude for gg ff via t- channel fermion exchange 
city amplitude for gg — > // via w-channel fermion exchange 
city amplitude for gg ff via s-channel gluon exchange 
city amplitude for gf — > /0 via t-channel scalar exchange 
city amplitude for gf — » /0 via t-channel scalar exchange 
city amplitude for gf — > /0 via s-channel fermion exchange 
city amplitude for gf —>■ f(f) via s-channel antifermion exchange 
city amplitude for gf —>■ fcf) via w-channel fermion exchange 
city amplitude for gf — > fcj) via w-channel fermion exchange 
city amplitude for gg ff via t- channel fermion exchange 
city amplitude for gg ff via w-channel fermion exchange 
city amplitude for gg f f via s-channel gluon exchange 
city amplitude for // — > // via t-channel gauge boson exchange 
city amplitude for // — > // via t-channel gauge boson exchange 
city amplitude for // ff via t-channel gauge boson exchange 
city amplitude for ff ^ ff via t-channel gauge boson exchange 
city amplitude for // ff via s-channel scalar exchange 
city amplitude for // ff via t-channel scalar exchange 
city amplitude for // — ^ // via w-channel scalar exchange 
city amplitude for // — // via s-channel scalar exchange 
city amplitude for // ff via t-channel scalar exchange 
city amplitude for // — ^ // via w-channel scalar exchange 
city amplitude for // —* ff via s-channel scalar exchange 
city amplitude for // ff via t-channel scalar exchange 
city amplitude for // ff via w-channel scalar exchange 
city amplitude for // // via s-channel scalar exchange 
city amplitude for // ff via t-channel scalar exchange 
city amplitude for // ff via w-channel scalar exchange 
city amplitude for // ff via s-channel scalar exchange 
city amplitude for // — > // via s-channel scalar exchange 
Table Subroutines. (Continues) 
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Name 


Description 


HWHSSl 


Gaugino-gaugino production matrix element 


HWHSS2 


Gaugino-gaugino production matrix element with polarizations 


HWHSSG 


Gaugino-gaugino/gaugino-sparton production 


HWHSSL 


Slepton pair production 


HWHSSP 


Combines MSSM subprocesses 


HWHSSQ 


SQCD 2 — > 2 hard subprocesses 


HWHSSS 


Identifies MSSM hard subprocess 


HWHVIJ 


Hard subprocess W/Z + jet production 


HWHV2J 


Master subroutine for all gauge boson + two jet processes 


HWHVVJ 


Dummy WW ,W Z ,Z Z production subroutine (see sect, 


HWHWEX 


Top production by W exchange 


HWHWPR 


Hard subprocess: W production 




Soft minimum-bias or underlying event 


HWMEVT 


Generates min bias or soft underlying event 


HWMLPS 


Generates longitudinal phase space 


HWMNBI 


Computes negative binomial probability 


HWMULT 


Chooses min bias charged multiplicity 


HWMWGT 


Calculates weight for minimum bias events 




Random number generators 


HWRAZM 


Randomly rotated azimuth 


HWREXP 


Random number: exponential distribution 


HWREXQ 


Random number: exp. dist. with cutoff 


HWREXT 


Random number: exponential transverse mass 


HWRGAU 


Random number: gaussian 


HWRGEN 


Random number generator (I'Ecuyer's method) 


HWRINT 


Random integer 


HWRLOG 


Random logical 


HWRPIP 


Random primary interaction point 


HWRPOW 


Random number: power distribution 


HWRUNG 


Random number: uniform + gaussian tails 


HWRUNI 


Random number: uniform 




Spacelike branching of incoming partons 


HWSBRN 


Generates spacelike parton branching 


HWSDGG 


Drees-Grassie gluon distribution in photon 


HWSDGQ 


Drees-Grassie quark distribution in photon 


HWSFBR 


Chooses a spacelike branching 


HWSFUN 


Hadron structure functions 


HWSGAM 


Gamma function (for structure functions) 


HWSGEN 


Generates x values for spacelike partons 



Table Subroutines. (Continues) 
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Name 


Description 


HWSGQQ 


Inserts g ^ qq part of gluon form factor 


HWSMRS 


Subroutine for MRST PDFs 


HWSSPC 


Replaces spacelike partons by spectators 


HWSSUD 


Sudakov form factor /structure function 


HWSTAB 


Interpolates in function table (for HWSSUD) 


HWSVAL 


Checks for valence parton 




Miscellaneous utilities 


HWUAEM 


Running electromagnetic coupling constant 


HWUAER 


Real part of photon self-energy 


HWUALF 


Two-loop QCD running coupling constant 


HWUANT 


Finds a particle's antiparticle 


HWUATS 


Replaces & with 


HWUBPR 


Prints branching data for last parton shower 


HWUBST 


Boost event record to/from hadron-hadron c.m.f. 


HWUCFF 


Coefficients for e^e~ and DIS cross sections 


HWUCI2 


Logarithmic integral Ci2 


HWUDAT 


Particle properties (N.B. BLOCK DATA) 


HWUDKL 


Generates decay vertex of unstable particle 


HWUDKS 


Converts decay modes into internal format 


HWUDPR 


Prints particle properties and decay modes 


HWUECM 


Centre-of-mass energy 


HWUEDT 


Insert or delete entries in the event record 


HWUEEC 


Computes coefficients for e+e~ cross section 


HWUEMV 


Moves entries within the event record 


HWUEPR 


Prints event data 


HWUFNE 


Finishes an event 


HWUGAU 


Adaptive gaussian integration 


HWUGUP 


Run termination for Les Houches Accord 


HWUIDT 


Translates particle identity codes 


HWUINC 


Initial parameter-dependent calculations 


HWUINE 


Initialises an event 


HWULB4 


Boost: rest frame to lab, no masses assumed 


HWULDO 


Lorentz 4- vector dot product 


HWULF4 


Boost: lab frame to rest, no masses assumed 


HWULI2 


Logarithmic integral Li2 (Spence function) 


HWULOB 


Lorentz transformation: rest frame lab 


HWULOF 


Lorentz transformation: lab rest frame 


HWULOR 


Multiplies by Lorentz matrix 


HWUMAS 


Puts mass in 5th component of vector 



Table P9| : Subroutines. (Continues) 
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Name 



Description 



HWUMBW 
HWUMPO 
HWUMPP 
HWUNST 
HWUPCM 
HWUPUP 
HWURAP 
HWURES 
HWUROB 
HWUROF 
HWUROT 
HWUSOR 
HWUSPR 
HWUSQR 
HWUSTA 
HWUTAB 
HWUTIM 



Generates mass (Breit-Wigner distribution) 
Spinor routine 
Spinor routine 

Converts integer to character 

Centre-of-mass momentum 

Prints contents of Les Houches common block 

Rapidity 

Computes/prints resonance data 

Rotation by inverse of matrix R 

Rotation by matrix R 

Computes rotation R from vector to z-axis 

Sorts an array in ascending order 

Print contents of spin correlations common block 

Square root with sign retention 

Makes a particle type stable 

Interpolates in a table 

Checks time remaining 



Vector manipulation 



HWVDIF 
HWVDOT 
HWVEQU 
HWVSCA 
HWVSUM 
HWVZRI 
HWVZRO 



Vector 
Vector 
Vector 
Vector 
Vector 
Vector 
Vector 



difference 
dot product 
equality 
times scalar 
sum 

zero (integer) 
zero 



HWWARN 



Issues warnings and deals with errors 



Table 2S: Subroutines. 



In addition there are the routines for generating the Schuler-Sjostrand parton 
distributions of the photon: 

SASAND SASBEH SASDIR SASGAM SASVMD 

and for putting the output of TAUOLA [88] into the /HEPEVT/ common block: 

FILHEP 

Finally, there are dummy versions of external routines, which should be deleted 
if the relevant packages are used: 

• CERN Library routine giving TRES= CPU time (in seconds) remaining (dummy 
returns TRES=10^°): 
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TIMEL 

• CERN PDFLIB structure function package: 

PDFSET STRUCTM 

• EURO DEC B decay package: 

EUDINI FRAGMT lEUPDG IPDGEU 

• CLEO B decay package: 

DECADD QQINIT QQLMAT 

• HERBVI baryon number violation package: 

HVCBVI HVHBVI 

• CIRCE beamstrahlung package (see below): 

CIRCEE CIRCES CIRCGG 

• TAD OLA T lepton decay package (see below): 

DEXAY INIETC INIMAS INIPHX INITDK PHOINI PHOTOS 

• Les Houches interface (see below): 

UPINIT UPEVNT 

9. Interfaces 
9.1 Les Houches 

Version 6.5 includes support for the interface between parton level generators and 
HERWIG using the Les Houches Accord as described in [89]. In general we have tried 
to code the interface in such a way that from the user point of view it behaves in the 
same way as that already included in PYTHIA [90]. 

The interface operates in the following way. If the IPROC code is set negative 
then HERWIG assumes that the user wants an external hard process using the Les 
Houches accord. 

If this option is used the initialisation will call the routine UPINIT to initialise the 
external hard process; this name is the same as that used by PYTHIA. This routine 
should set the values of the run parameters in the Les Houches run common block. 
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After the initialisation during the event generation phase the routine UPEVNT 
(again this name is the same as that used by PYTHIA) is called. This routine should 
fill the event common block as described in [89] . 

Dummy copies of both these routines are supphed with HERWIG and should be 
deleted and replaced if you are using this option. Due to the internal structure of 
HERWIG two new parameters are needed to control the interface, in addition to those 
in the Les Houches common block. The logical input variable LHSOFT [ . TRUE] controls 
the generation of the soft underlying event; the default is to generate a soft underlying 
event and so LHSOFT must be set .FALSE, for lepton-lepton processes. The second 
variable LHGLSF [.FALSE.] controls the treatment of colour self-connected gluons, 
which may occur with some methods of colour decomposition. The default is to kill 
events with self-connected gluons whereas if LHGLSF = . TRUE . an information-only 
warning is issued instead. 

We have not included full support for the interface. In particular 

• It is assumed that all events have two particles of status IDUP=-1; 

• We do not support the status codes -2 and ~9; 

• The treatment of the code IDUP=3 is the same as IDUP=2, i.e. all intermediate 
resonances have their masses preserved. 

These restrictions are the same as those imposed by PYTHIA. If you have any prob- 
lems with the interface or need the options we have not yet supported please let us 
know. 

It should also be noted that while the interface has been tested with Standard 
Model processes, for example all the processes in ALPGEN [91], it is less well tested 
for SUSY processes. 

9.2 CIRCE 

Simulation of beamstrahlung is now available, via an interface to the CIRCE pro- 
gram [92]. It is implemented as a modification of the standard bremsstrahlung im- 
plementation, so is only available for processes for which that is available. The 
produced radiation is treated in the coUinear limit, i.e. COLISR and the option WWA 
in routine HWEGAM are forcibly set . TRUE . . 

For both types of distribution function, /e/e and /^/e, we use a simplifying ap- 
proximation: where the correct result should use the convolution of beamstrahlung 
and bremsstrahlung, we actually use the sum of the two. For the photon distribution 
this is quite reasonable, but for the electron it is perhaps more questionable. It boils 
down to the replacement: 



e/el-^) / /e/e,brem ( -2- j /e/e,beam ( ) ^ /e/e,brem(3^) ~l~ /e/e,beam(3^ 






(9.1) 
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This is good to the extent that both distributions are dominated by the x —>■ 1 
region. A small utility program can be obtained from the HERWIG web page to test 
this approximation. For example, for TESLA running at 500 GeV, the integral over 
fe/e{x) for large x, relevant for threshold shapes, is accurate to 1% for x > 0.95 and 
to 10% for X > 0.99, and the value of fe/eH-^Y), relevant for radiative return to the 
Z, is accurate to 2%. 

In order for the distribution to be probabilistic, the coefficient of the delta func- 
tion has to remain positive. If it is not, HERWIG will terminate. This can be 
prevented by adjusting the ZMXISR parameter: the default value, ZMXISR=0. 999999 
is too large and HERWIG will fail, but reducing it to ZMXISR=0. 99999 works fine 
for all the parameter sets currently available. It is straightforward to check that 
ZMXISR values in this range do not have an infiuence on any physical observables, by 
rerunning with it further reduced. 

The interface is controlled by five new variables, which are given in the following 
table together with their default values: 



CIRCOP 





CIRC AC 


2 


CIRCVR 


7 


CIRCRV 


9999 12 31 


CIRCCH 






CIRCOP is the main control option: CIRCOP=0 switches off beamstrahlung and uses 
standard HERWIG; CIRC0P=1 switches to collinear kinematics, but still leaves beam- 
strahlung switched off; CIRCGP=2 uses only beamstrahlung; and CIRC0P=3 uses both 
beamstrahlung and bremsstrahlung. CIRCOP=0 and CIRCGP=3 should therefore be 
regarded as 'off' and 'on', respectively, with the other two options mainly for cross- 
checking purposes. The variables CIRCAC, CIRCVR, CIRCRV and CIRCCH are simply 
passed to CIRCE as its input variables acc, ver, rev and chat, as described in its 
documentation. The default values correspond to the most up-to-date revision of 
version 7 of the TESLA parametrization, with minimal output. 
CIRCE can be obtained from 

http : //heplix . ikp .physik . tu-darmstadt . de/nlc/beam.html 
9.3 TAUOLA 

An interface to the TAUOLA decay package [88] has been added. To use this interface, 
the dummy subroutines DEXAY, INIETC, INIMAS, INIPHX, INITDK, PHOINI, PHOTOS 
must be deleted. You should then link to both TAUOLA and PHOTOS. The easiest 
way to do this is to obtain the TAUOLA/ PHOTOS versioning system from 

http : //wasm . home . cern . ch/ wasm/test . html 
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and produce a version with the correct size of the HEPEVT common block. This system 
will provide a JETSET demo. If you link HERWIG with the code produced by this 
apart from the JETSET interface and main program everything should work. 

This interface uses the information from the spin density matrices to select the 
helicity of the decaying r leptons if the spin correlation algorithm is being used; 
otherwise the helicity of the decaying r is averaged over. 

9.4 ISAWIG 

ISAWIG is an interface which allows SUSY spectra and decay tables generated by 
ISAJET [39] to be read into HERWIG. The web address is 

http : //www.hep .phy . cam. ac .uk/~richardn/HERWIG/ISAWIG/ 

To coincide with the release of HERWIG 6.5 we have produced a new version 
of ISAWIG. Since the original version of ISAWIG there have been a number of new 
versions of ISAJET, often with changes in the sizes of the common blocks from which 
we extract the information we need. This has necessitated periodical updating of 
the code to run with the most recent version of ISAJET, with the result that the 
code could no longer be run with older versions. However, the Snowmass points and 
slopes [93] are defined with ISAJET version 7.58 and so we can no longer continue 
in this way and still be able to generate these points. Therefore, starting with the 
new ISAWIG version 1.2, we are using C preprocessing so that users can define the 
version of ISAJET they are using at compile time. When compiling the main ISAWIG 
code and the modified SUGRUN and SSRUN programs the following compiler options 
should be specified: 

• -DISAJET758 to use ISAJET7.58, 

• -DISAJET763 to use ISAJET7.63. 

The default at the moment is to compile code to run with ISAJET7.64. In the future 
this will change so that the most recent version of ISAJET becomes the default but 
we will continue to support the older versions. 

9.5 HDECAY 

There has also been a new release of the HDECAY package [94], version 3. In order to 
support this version, and the previous version 2.0, we are also using C preprocessing 
to control the version of HDECAY used. In order to achieve this the HDECAY interface 
code has been merged with the main ISAWIG program. The following options should 
be used when compiling ISAWIG if you are using HDECAY: 

• -DHDECAY2 if you are using version 2 of HDECAY 

• -DHDECAY3 if you are using version 3 of HDECAY. 
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As before the default is to compile a dummy routine. If you are using HDECAY you 
must link with a version of the HDECAY code which has the main HDECAY program 
removed. 



9.6 MC@NLO 

The program MCONLO [95] generates events with 

• exclusive rates and distributions accurate to next-to-leading order when ex- 
panded in as] 

• multiple soft /collinear parton emission generated by HERWIG parton shower- 
ing; 

• hadronization according to the HERWIG cluster model. 

At present only the processes of gauge boson pair production in hadron [h = p ot 
p) collisions are implemented in MCONLO, making use of the HERWIG process codes 
IPRGC=2850-2880 (see table |). 

To use these codes, one must delete the dummy HERWIG subroutine HWHVVJ and 
replace it by the one in the MC@NLO package. This will then read and process an 
input file of parton configurations prepared by MCONLO. For further details see the 
MCONLO web page: 

http : //www . hep . phy . cam . ac . uk/ theory/webber/MCatNLO/ 

Implementation of other processes, such as heavy quark production, in MCONLO 
is in progress [96]. These future development will make use of the Les Houches 
interface discussed in sect. |9J] , rather than special process codes. 

9.7 JIMMY 

Contrary to our original plans, the JIMMY generator [49] for multiple parton-parton 
interactions in hadron-hadron, hadron-photon and photon-photon collisions is still 
not integrated into HERWIG, and is still a separate add-on package. It can be ob- 
tained from the JIMMY web page: 

http : // j etweb . hep . ucl . ac . uk/ JIMMY/ index . html 

where references to the relevant papers can also be found. 



9.8 HERBVI 

An interface is provided to the HERBVI package for electroweak baryon number vio- 
lation (p), and other multi-W^^ production processes [44]. For full details, see 
http : //www ■ hep . phy . cam . ac . uk/| r^richardn/HERWIG/HERBVI/ 
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